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AN ANALYSIS OF A CHARRING ABLATION 


THERMAL PROTECTION SYSTEM 

By Donald M. Curry 
Manned Spacecraft Center 


SUMMARY 


An analytical model is presented for predicting the transient one- 
dimensional thermal performance of a charring-ablator heat-protection system 
when exposed to a hyperthermal environment. The heat-protection system is 
considered to consist of an ablation material and backup structure. The 
ablating material is further considered to consist of three distinct regions 
or zones: char, reacting, and virgin material. 

A FORTRAN IV digital computer program (STAB II) utilizing an implicit 
finite difference formulation has been written for the IBM r J09b/bO computer 
system. The program considers one ablating material and a maximum of 12 back- 
up materials with conduction or radiation and/or convection allowed between 
materials. Thermal properties of all materials are temperature dependent, 
with the. properties of the charring material also being state dependent. 

The governing differential equations and their implicit finite difference 
formulation are presented. The program input and output are described in detail. 
The FORTRAN program statements and nomenclature are presented. Also, the 
theoretical and experimental results are compared. 


INTRODUCTION 


The analysis and design of thermal, protection systems for entry into an 
atmospheric environment have resulted in a voluminous amount of literature on 
the general subject of ablation. (See refs. 1 and 2 for a survey of information 
on ablation.) The ablation materials may generally be classified into three 
categories: subliming, melting and vaporizing, and charring. The charring 

ablator normally provides the most efficient thermal-protection shield for the 
major portion of a manned entry vehicle. This report' describes a method for 
predicting the thermal response of a typical charring-ablation material. The 
response of a charring material to a hyperthermal environment is 'extremely 
complex, and the mathematical model presented to analyze the transient behavior 
of the material contains simplifying assumptions and approximations necessary 
to afford even a numerical solution. 



The equations derived in this analysis have "been programmed in FORTRAN 
IV for an IBM r J0 ( jb/k-0 computer system. The numerical formulation of this 
digital program, designated STAB II, is such that an implicit solution is 
obtained. The thermal response of a typical charring material as predicted by 
STAB II Is compared with arc tunnel results. 

A sample problem is presented in appendix A. Program usage instructions, 
including definitions of the input terminology, are presented in appendix B. 
Appendixes C and D are the program FORTRAN IV statements and definitions of 
the program terminology. A general flow chart of the program is presented in 
appendix E. 


SYMBOLS 


A 

c 


P 


E 

F 

F 


env 


H 


d 


H v 

h 300 

h 

h 

env 


collision frequency 
specific heat 

activation energy 
exterior view factor 

view factor— emissivity product to cabin environment 

heat of virgin material degradation 

total enthalpy 

wall enthalpy 

enthalpy of air at 300° K 

film coefficient between backup materials 

film coefficient between last backup material and cabin 
environment 


k thermal conductivity 

fa mass loss rate of char material 

c 

fa gas ablation rate 

g 

NP number of nodes in ablation material 
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c blow 


order of reaction 

net heat rate into front surface 

hot wall convective heat flux with blowing 

heat flux due to combustion 

cold wall convective heat flux without blowing 

radiation heat flux 

universal gas constant 
surface recession depth 
surface recession rate 

temperature of node at beginning of time step 

cabin environment temperature 

temperature of node at end of time step 

radiation heat sink temperature 

thickness of ablation material 

distance from surface to any point 

heat of combustion per unit weight of char 

thickness of a node 

time step (9' - 0) 

emissivity of material 

transpiration cooling efficiency 

initial time 

final time 

transform for the ablation material 
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p density 

a Stefan- Boltzmann c ons tant 

\|f blocking effectiveness function 

Subscripts: 

c charred state 

i node number 

j material number 

v virgin state 


PROGRAM DESCRIPTION 


The following general requirements were established before writing a 
digital computer program to analyze a charring ablation system: 

(1) Stability of the equations for all applications. 

(2) Machine running time short enough to make use of the program 
economically feasible (a minimum of turn around time per problem). 

(3) A minimum of input per problem. 

(4) A wide variety of boundary conditions for application to both 
trajectory data and ground or flight test data analysis. 

STAB II has been formulated in FORTRAN IV to analyze the transient 
thermal performance of a charring ablator heat protection system. The program 
considers one ablating material and up to 12 different backup materials with 
or without air gaps. Pure conduction or radiation and/or convection between 
backup materials is allowed. The ablation material may be divided into a 
maximum of 50 nodes, and each backup material may be subdivided into a maximum 
of 10 nodes. The thermal properties of the materials are in tabular form and 
are temperature dependent. The ablation material is also dependent upon its 
state, that is, fully charred, partially charred, et cetera. 

The following surface boundary condition options are provided: 

(l) Cold- wall convective and radiative heat flux tables as a function of 
time. These components are specified separately, since mass transfer at the 
surface blocks part of the convective heating but, in general, has no effect on 
the radiant heating. 
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(2) Surface temperature as a function of time. 

(3) Surface recession as a, function of temperature or time. Surface 
recession as a function of temperature and pressure is also available. 

Heat loss to the interior environment for the last node of the backup 
structure can be specified by two methods: 

(1) Conduction into the node and radiation and/or convection loss to the 
interior environment. 

(2) Conduction into the node and adiabatic wall. 

The STAB II numerical formulation of the equations describing the response 
of the heat shield is such that an implicit solution has been obtained. It is 
well known that numerical solutions of partial differential equations are 
subject to several different types of errors. The first of these is the trun- 
cation error, due to the use of a finite subdivision. This error may be reduced 
by simply choosing a smaller subdivision, AX. The exact values are approached 
more and more closely as AX decreases. The second kind of error is the nu- 
merical, or roundoff error. The way in -which this numerical error grows or 
decays with time determines the stability of the difference equations. 

To il lustrate the differences in the explicit and implicit equation form, 
consider a nonablating homogeneous solid. The one- dimensional Fourier conduc- 
tion equation, neglecting any heat generation terms, is 



The finite difference form of equation (l) written in the conventional forward 

. th 


time step or explicit form for the i 

w 


node is 


cf 

. o 


Sfcfr (T tl - * t ) (T t - T 1+1 ) ax (t; - Tj 

fV- . AX . AX " AX . AX pc p A0 


2k 


i-1 


AX 

2k. 


AX 

2k. 


( 2 ) 


2k 


i+1 


where the prime superscript denotes values at the end of the time step 

A0 = e l - e 

For explicit conduction solutions, the following stability criterion has 
been established: 

PC P (AX)f > 2 

k A0 
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which places an upper limit on the time step A0 for a fixed truncation error. 
This criterion can require a prohibitive amount of machine time. 

Liebmann (ref. 3 ) advocated a solution of the equation which does not 
require this stability criterion. The finite difference equations are written 
in a backward time step form which affords an implicit solution. 

The implicit (backward time step) difference form of equation (l) for the 
i^* 1 node is: 



AX AX_ " AX_ 

2k. . 2k. 2k. 

l-l l l 



AX 


2k 


1+1 


Pc 


AX 



( 3 ) 


Equation ( 3 ) uses the temperature differences at the end of the finite time 
interval instead of the beginning, as in the explicit method of equation (2). 
The only known temperature in equation ( 3 ) is T^, but there are corresponding 

equations for each point in the system, and all are solved simultaneously to 
yield the temperature at each node. 

Collecting all unknown temperatures on the left side of the equation and 
the known temperature on the right side, equation (3) becomes 


/ 1 


\ T U-( 

£X . 

AX 

K-i 

•H 

CVI 

1 l 

c AX\ 

p i L 

1 j- / 

f i 


AX AX_ AX_ AX 

2k. . 2k. 2k. 2k. 

l-l l l l+l 


'p.c AX\ 


AX AX 
2k. 2k .. 

l 1+] 


\ T t 

'i+1 


1 p i 


A^yk 


Equation ( 4 ) is of the form 


AT’ + BT' + CT' = D 
i-1 i i+1 


STAB II generates such an equation for each node in the system. 




( 5 ) 


Since radiation is an important mode of heat transfer in charring ablative 
systems, a problem is encountered in any equation containing a radiation term. 
The radiation heat flux, written in a backward difference form is: 


^■rad 


Fea (t^ 4 



(6) 


*J 

I t? 
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This term cannot be used in an implicit solution since the unknown temperature 
T_J is to be the 4th power. The 4th power unknown can be eliminated by the 

following linearizations: 


where 




+ AT 





let 


AT = T* - T 
i i 


Z 


AT 

T. 

1 


and rewrite equation (7) as 

( n) 1 * - ( + 

If Z has an absolute value near zero , the following is true 

(l + zb £ 1 + kz 

Now substituting equation ( 9 ) into equation (8) 


( t ;) 4 3 ( T i) 4 ( i+4z > -( T i) 4 ( i+4 f) 





( 7 ) 


(8) 

(9) 


(10) 


Equation (10) is a linearized approximation of equation (j) in which the unknown 
temperature is only to the first power. The assumption in equation (10) is that 
AT/!, has an absolute value near zero. Figure 1 is a plot of the error obtained 
1 l\. ' 

when (l + 4Z) is substituted for (l + Z) . For most ablation problems in which 

the surface temperature is high and the radiation losses are significant , the 
value of AT/T can easily be controlled to values of less than ±0.1. 


Therefore, 


equation (6) 


Srad 


can now be written 
Fea^T^T^ - - 



( 11 ) 
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Using the linearized approximation for the radiation terms , the resulting system 
of implicit difference equations constitute a tridiagonal matrix of the follow- 
ing form: 


B 1 T 1 + C 1 T 2 


= D, 


A 2 T 1 + B 2 T 2 + C 2 T 5 


= D, 


A 3 T 2 + B 3 T 3 + C 3 T 4 


= D, 


^N-l + B N T N D N 


Gauss T elimination method, discussed in reference 4, is applied to solve 
the system of equations. This method affords a fast and accurate solution for 
matrices containing a dominant diagonal. The solution of_thi s matrix gives the 
t emperature of each f or^tHe^pext future. time step- The 

entire process is repeated for each time step throughout the run, giving a time 
history of the temperature at each node. 


Using this method of solution, residual errors in the temperature compu- 
tations at the beginning of the time step are distributed throughout the entire 
system of nodal equations and tend to cancel out rapidly. The principal advan- 
tage in using the implicit method is a set of equations that are mathematically 
stable in time and distance. Therefore, the magnitude of the time step is not 
limited by a convergence criterion. However, care must be taken in selecting 
the magnitude of the time step in order to minimize truncation errors when the 
second derivative of temperature with respect to time is large. A similar 
approach is used to minimize truncation errors in distance by choosing small 
node dimensions in locations where large second derivatives of temperature with 
respect to distance are expected. 

In the case of a char-forming ablative heat shield where approximately 
80 percent of the heat is reradiated, instability can arise in taking large 
time steps. The temperature of the surface node can start oscillating on suc- 
cessive time steps when a balance between the radiation source and the heat sink 
has been achieved. Therefore, in ablation problems in which the surface node 
loses a large percentage of heat by radiation, oscillations of the node can be 
damped out by taking small time steps during* conditions of high heat flux and 
near radiation equilibrium temperatures. 


ANALYSIS 


Figure 2 is a schematic of the thermal protection system to‘ be analyzed. 

A receding surface has been assumed with the formation of a residual char layer 
and reaction zone. The thermal protection system is composed of one charring 
material and a maximum of 12 different backup materials with or without air gaps. 
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The analysis is such that the entire system may be composed of noncharring 
materials. The thermal properties of all materials are temperature dependent; 
also, the charring material properties are state dependent (fully or partially 
charred) . 

The response of charring ablation heat shields to a hyperthermal environ- 
ment is extremely complex, and simplifying assumptions and approximations are 
necessary to afford a numerical solution. The following assumptions and approx- 
imations are utilized in the equations developed in this report: 

(1) The material decomposes from the virgin state to a porous char layer 
in the reaction zone. 

(2) The reaction zone can be defined by an upper and lower temperature 
limit. 

(3) The gas generated within the reaction zone is assumed to pass out of 
the structure with no pressure loss. No gas accumulation within a node is 
allowed. 

(k) Local thermal equilibrium is maintained between the gas and porous 
char matrix. 

(5) The gas undergoes no further chemical reaction within the residual 
material after having been formed. 


Derivation of Equations 

The equations are derived for a moving boundary coordinate system, where 
the front face is the moving surface (ref. 5). With this system, the ablating 
material is divided into a fixed number of nodes of thickness AX which depends 
on the instantaneous location of the front face. The surface recession is 
handled in a continuous manner, eliminating the need of throwing away or lumping 
off of nodes. 

The physical model for the front surface, including all heating terms, is 
shown as follows: 
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AX 

2 . 0 


UU; 

4 - C. ^ 


Q. 

m 


m c T. 
§1 Pi i 

Vp/A 



tfv I 


The energy equation at the front char surface is 

dT 


& (i ^ p i c Pl Tl , 


= | AX Pl c pi ^ + | 


= Q. + m c T^ + P„S^c T_ - m c T_ 

m g 2 p 2 2 T 2 p 2 2 g x p x 1 


• 1 /at\ 

-Pc ST - k I— I 
P l C Pl 1 X 1 1-2(AX / 


(12) 


where 

^in ~ *e, blow + ^rad + ^eomb “ Fecr (^1 ^=° ) 

and 

d(AX) _d_ f VL - S \ S 

d§ de \KP - 1 J " NP - 1 

where S is the linear surface recession rate and UP is the total number 
of nodes in the ablation material off thickness ^VL."] 
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xn 


finite difference form 

Sp,c TJ - 
1 P 1 1 

( T i - T a) 

AX , AX 
2k 1 2k 2 

n ' AX( 

— * f) /I r \ 

T* - T \ 
1 1 1 

1 Pi 2 \ 

A0 ) 

- Kvi 

{ s ) 
VHP - 1/ 


(12a) 


Then, rearranging and collecting terms yield 


me + Sp c + p c — - + -—pc ( - \ 

g p, M 1 p M 1 p 2A0 + AX_ AX_ 2 P l C p I MP - 1 I 

- ~ ^ - 1 - Pk Pk x ' / 


2k l 


me + 


vp 2 ■ « ♦ + P2 \ s 


2k i 


2k, 


sfo - 

l NP - 1.0 J 


T* 


f r ^’ c 


" P l C Pl 2A9 T 1 ~ Q in. 


(12b) 
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The physical model for interior points in the mature char zone, including 
all heating terms, is shown in the following sketch: 



The energy equation for interior points in the char matrix is 



AX. p.c 
i i P i 



= AXp.c 
1 p 


dT. 

1 


d£ 


P.c t! 
i p ± i 



-0 


= m c T. + k 
g i+l p i+l 1+1 


/at \ 

1-1, i ^AX / 


. (W? - i - 

+ Pi+lC p i+1 S V wp - i / T i+i 


- k. - m c T| - p.c S (— 1 - g ) t' (13) 

i,i+l \ AX y % Pi 1 ip. \ HP - 1 / i 
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Putting equation (13 ) in an implicit finite difference form yields 



( 13 a) 


In the mature char zone, no internal gaseous ablation products are assumed 
to form. The reaction zone is the source for the formation of the internal 
gaseous products. Therefore, in equations (12) and (15^ = m 

S i S i+1. 

The physical model for nodes in the reaction zone is identical to sche- 
matic shown for the interior nodes in the char except for considering the energy 
absorbed in formation of the gaseous ablation products. The heat balance 
equation for a node in the reaction zone is 


-tx(/SK p.c T.\ = AXp . c ry 1 - P.c T l ( — - ) - An - m 

dBI 1 p ± i) 1 p ± d0 i p ± x \NP - 1/ ( g ± i 


v) 


= m c T. , n + k 
S i+1 p i + l 1+1 


L T i+l + k i-l,i (ax) + P i+l c p 1+1 ^ \ UP - 1 / T i+1 

- k i,i+l (if) ■ 4, g i °P 1 T i ■ P i°P 1 S V HP - 1 / T i 
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I I III 


■■nil 


III II II II II 



Rearranging, 



The physical model for the interface between the reaction zone and virgin 
material is illustrated as follows: 


— AX ► 
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The heat balance equation for this node is 
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The physical model for an interior node in the virgin material is 



- k i-i,i(S) + p i+i c P s 


TCP - 1 



Rearranging, 



The physical model for the last node in the ablation material and first node 
in the backup structure is 



For this interface, T. . = T. . _ 

i, 3 i, J+l* 
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The heat balance equation for this node is 


d_ 

ae 


/AX. AX \ 


( AX.P. .C + AX..-C P. . A , . 

0 1,0 Pi J 0+1 Pi „1+1 l,0 +1 ) dT i 1 / s \ . 

J 2 jsr-s{m> —) \ , p i,j T i 

± 9 

0 / AT \ „ a I 2 1 m . , /at\ 

' k i-l,J (AX ) " C p. P i,j S \NP - 1/ i “ i, j+1 (ax) 

n* ^ f ± 9 J n4-1 n4-l ' 


(17) 




Rearranging yields 


T\ 


AX. AX. i-1 

J + <L 


2k. . ‘ 2k. . 

i-lj. J i, J 


i+1, j+1 


+ + L 


AX. AX. T AX. . AX. n 

J j. . J_ ,1+1 , ,1+1 


2k. _ . 2k. . , 

1 1,J 2k i,j ^ J+1 2k i+l,j+l 


fAX.P. .c + AX. t ,c p. . n , 

j p, ; . J+1 p lf - +1 i,J+ U 

2A0 


T ' +■ 

i f AX. 


J+ 1 + ~j+l 


rp 1 

AX.., Vi+1 


2k. 2k. 

1, J+1 1+1 > J+ly 


t 


AX .c p . . + AX . c 

J P, 1 11 > J J+1 P, 

±2 v x 

2A9 




T. 


(17a) 


The backup struct ure may contain up to a maximum of 12 different materials 
with or without air gaps between materials. Therefore, conduction or radiation 
and/or convection between materials is allowed. The heat balance equations for 
the various modes of heat transfer in the backup structure are presented in 
the following equations: 


(l) Interior node material: 


T« 

l 1-1, J 

- t; .] 
1,0/ 

T. 1 . 
1 i,J 

- TJ .] AX. 

1+1,0/. - DC J 

AX. 

0 

aX. 

+ J 

AX. 

J 

AX P i,0 p A0 

a. J 1, J 

2k. . . 

i-l ,0 

2k. . 
i, J 

2k. . 
i, J 

2k. , . 
i+l, 0 


Rearranging, equation (l8) "becomes 


t! 


AX. 


AX. i-l,j 


AX. 


AX. 




2k. . 

1,0, 


AX. AX. 

a_ * a. 


v2k 




2k. . 2k. . • 2k. . , 

i ,0 i ,0 1 + 1,0 


(18) 


+ P, ,c 


AX.\ 

a. \ rp I + 

v lfi A9 1 i,J 


AX. 


t : , 


ax. ri+i,j 


2k. . 2k. ^ . 

1,J 1+1,0, 


AX, 


= -Pi ,c 


.« tr t - • (l8a) 

,J P, , A6 i,0 
- 1 ? J 


(2) First and last nodes of two interior materials with no gap: 


m T m I 

1 -1>3 ili 

AX. AX. 

0 + 1_ 

2k. . . 2k. . 2k 

i-l, J 1,0 


rp I m 1 

i,o+i " i+i,o+i 


U ~ _L 

Olr 

' 1 , 0+1 1 + 1 , 0+1 


p. .c AX. + P. . _c AX 

i,J P lti J i,J + lp. t . +1 0+1 


2 A0 / V T i,j “ T i, j 


(19) 


For this case. T. . = T. . , _ 

’ 1, 0 1, j+l 
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Rearranging, equation (19) becomes 

^r\ T i-i,j 

J JL lL 


2k - 2k , 
i-l,J x,j 


AX. AX. AX AX 

J_ + sL_ J+l- , 3+1 

2 k. . . 2k. . 2k. ... 2k. _ ... 

i-l, J i,J 1,0+1 i+l,j+l 


/p. .e AX. + p. c AX.,-\ 

+ 1 1,3 P i,.i J X * J 1 P i.,i+1 3+1 j 

' 2 A9 ' 


T! . 


**,1+1 , 5l+l ^ i+1 ^ +1 

2k i,J+l 2k i+l, j+l 


/ p. .c AX. + c p. AX \ 

J+1 L 

V 2 AS / i,j 


(3) First node of interior material -with an air gap between materials: 


(19a) 


h 


j( T i-l,j “ T I,j+l) + ^j- + -i- - ij V*' 1 ’* " T ^J +1 ) 


/ 

m * m* 

1..1+1 ~ i+l,. i+l 
AX AX. 

J±i_ + J+ l . . 

2 k. ... 2 k.,_ ... . 

^ i, J+l i+l, J+l / 


p. . - c AX 

1 ’ J+1 Jtl 

2 AG 


m ' m 

l i,J+l " i,j+l 


(20) 


Equation (20) may be linearized by using the approximation 

k x k 

t T ~ T - yjy 


as discussed in the Program Description section. 


Therefore, rearranging and linearizing, equation (20 ) becomes 



P. . , c 
.1, J+l P 


AX 


i> j+l 


j+l 


AX 

. 1+1 

2k. 

i, J+l 


AX 


j+l 


2 A0 


T t 

1 i, j+l 


pv 

i+1, j+l 


1 

AX . ' AX 

■1+1 + J+l 


2k. ... 

i, j+l 


2k 


i+1, j+iy 


T ' 

1+1, j+l 


P. 

i, J+l P 


j+l 


AX. n 
J+l 


2 A0 


T 

1, j+l 



(4) Last node of an interior material with an air gap between materials: 


t! - t! . ] / \ 

i~ 1 1 j 1? j 1 _ v, I mi _ T' l 

AX. AX. j\±,j l,j+ly 


AX. 

J 

2k. , . 2k. . 

l-l, j i,j 



21 



Rearranging and linearizing , equation (21 ) becomes 




( 5 ) Final node in backup structure: 
(a) Adiabatic surface — 



Rearranging, equation (22) becomes 



(21a) 


( 22 ) 



(22a) 
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(b) Radiation and/or convection loss to cabin environment — 



Discussion of Assumptions 

A brief discussion of several assumptions and approximations made in 
deriving the heat balance equations is now presented. 

As shown in the Derivation of Equations section* transient heat con- 
duction* thermal degradation, and the flow of the gaseous products from the 
reaction zone are the internal thermal transport phenomena of interest. 

Several methods are available in the treatment of the thermal decomposition 
process* and they differ primarily in whether the chemical decomposition 
occurs in a single plane at a fixed temperature or whether a spacially contin- 
uous decomposition in depth is assumed. This analysis assumes that the decom- 
position from the virgin to . the char state occurs in a reaction zoiieT'that ..i.s 
def ined by known temperature^ limits . These temperature limits are determined 
from thermogravimetric test data' for* the particular material being investigated. 
Figure 3 is a thermogravimetric curve for typical charring ablation material. 
From this curve* the rate of pyrolysis m is calculated by knowing the 
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temperature change of a particular node with time, that is. 


p i - 


p i - 


■aS" 


( 24 ) 


and 


m 




HP 

= E 

i 


P.AX. 


(25) 


This method of computing the gas-generation rates and local instantaneous 
density may he subject to error since the thermogravimetric curve of a 
material is influenced by temperature rise rate (deg/sec), and the reaction 
zone may shift up and down the temperature scale. This error can be elim- 
inated by the use of an Arrhenius expression of the form 




JL 

“rt 


( 26 ) 


^ The method now being used in STAB II (equation ( 25 )) t<^ calculate the pyrol ysis^ 
N^rate is being investigated to determine its validity^ The final formulation 
of the pyrolysis rate law must rest heavily on the experimental rate data for 
the material under investigation. The use of simple expressions such as 
equations (24) and ( 25 ) may be entirely adequate, depending upon activation 
energy for the decomposition process and order of reaction. 


The aerodynamic heating input in the analysis consists of convective 
and radiative components treated separately. This distinction is necessary 
since the convective heating can be significantly reduced as a result of the 
injection of the ablation gases into the boundary layer, with generally no 
effect on radiant heating. Reduction in the convective heating rate can be 
approximated by the following expression (ref. 6): 


^ block - 


H, 


T 


Therefore, 


H 


vj 


( 27 ) 


q , _ = q. 

c , Plow cw 


F T“ 


H 

w 


It “ 


loo, 


- q 


"block 


(28) 
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However, equation (28) is unsatisfactory for high blowing rates, since q block 

can become greater than q . An experimental curve of blocking effectiveness 
• cw 


q^ blow m 

i|r = — — as a function of the mass transfer parameter | - can be em- 


cw 


ployed to determine the heating reduction at high blowing rates. Both methods 
have been employed in the STAB II analysis. Equation (28) is presently in use. 
However, no satisfactory method for accurately predicting the convective heat 
blockage has been determined. 


Another source of heating is the combustion of the ablation products in the 
boundary layer. Reference 7 presents an analysis of the oxidation of a carbon 
surface and the resulting combustive heating. The heating due to combustion as 
derived in reference 7 is 


^■comb 


m AH 
c c 


(29) 


where AH is the heat of combustion per -unit weight of char, 
c 

The thermal properties of the ablation material are both temperature and 
state dependent (fully or partially charred). Figure 4 is an illustration of 
the variation of these properties with temperature and state. The thermal prop- 
erties are assumed to vary as follows: 


( 1 ) Char zone (t. £ 

k = f (temp) 


c = f (temp) 


p = constant 
c 


(2) Reaction zone (^51 


< T> < t 


char 


) 


p = f (temp) = P v + (P v - P c ) 


T i ~ T abl 
1 tibl “ lhar 


k = f(p) 


c p = f(p) 


= k + (k - k 

C \ -V c 


= c + /c - c 


He - p 
{\ ' \) 


V ’ Cj 

p_. - p. 
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(3) Virgin zone T abl 

p = constant 
v 

\ = f (temp) 

c = f (temp) 

P 

v 

The calculation of char removal, due to chemical, thermal , or mechanical 
mechanism or a combination of these mechanisms, has been examined by a multitude 
of investigators and numerous correlations exist, depending on the specific 
material involved. 

To provide a maximum degree of flexibility for analyzing both ground and 
flight test data and synthesizing trajectories, the following provisions for 
char removal (surface movement) are provided: 

(1) Removal of char as a function of surface temperature. 

(2) Removal of char at a rate which is a function of time. 

As the char is removed, the surface moves with respect to a coordinate fixed 
in the material. The distance between the initial surface location and the 
char surface is 


S 



ANALYSIS VERIFICATION 


As discussed in the previous sections, approximations and assumptions 
were made in the analytical model to afford a quick and accurate solution 
in predicting the thermal response of a charring heat shield. These 
simplifying assumptions and approximations are expected to introduce only 
minor errors; however, the validity of the analyses and resultant accuracy 
can be judged only by a comparison with exact theoretical solutions and 
experimental data. Three examples have been selected, and a comparison 
of the STAB II results with the theoretical and test data is discussed in the 
following paragraphs. 

An elementary transient heating example was chosen to demonstrate the 
accuracy and numerical stability of the STAB II program. A steel slab 
6 inchgs thick was selected and assumed to be at uniform initial temperature 
of 460 R (0° F). The thermal properties were considered constant. The 
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front surface was subjected to a heating rate of 72 Btu/ft -sec, and an 
adiabatic back surface was assumed. Figure 5 shows a comparison of the 
STAB II calculated in-depth temperatures as a function of time with the 
exact solution taken from reference 8. 

To demonstrate the STAB II solution with a moving boundary, a slab with 
constant properties, uniform initial temperature, front surface moving with 
a constant velocity, and constant surface temperature was chosen. The exact 
solution for a semi-infinite slab with these boundary and initial conditions 
is presented in reference 9- Figure 6 presents a comparison of the STAB II 
temperature response with the exact solutions. As can be seen from this 
figure, the two solutions are not in agreement for approximately the first 
50 to 60 seconds of the transient. This disagreement is the result of the 
quasi-steady state assumption made in the exact solution analysis. 



A calculation was made to estimate the induction time (time at which 


8t 

6e 


is a good assumption) and found to be approximately 60 seconds, which is 
agreement with the STAB II results. 


= 0 
in 


Finally, to verify the fully charring ablation model, an example of a 
typical charring material was chosen. (See the sample problem in appendix A.) 
The charring ablation material is initially 1.6 inches thick with an adiabatic 

2 

back surface and a constant heat flux of 95 Btu/ft -sec applied to the front 
surface. The surface is assumed to recede at a constant velocity of 

3.05 x 10~^in./sec. Figure 7 presents a comparison of the in-depth tempera- 
tures with actual test results obtained in an arc tunnel. The results are in 
good agreement, with the largest deviations between calculated and measured 
values occurring for the thermocouple located at a depth of 1.0 inch. The 
disagreement could be attributed to several possible errors: thermal property 

values, incorrect location of thermocouples, et cetera. The effect of varying 
the thermal properties (thermal conductivity, specific heat, et cetera) is 
presently being investigated. 


Tables I and II present the input and output data used for this example. 
Figures 8, 9? and 10 are the resulting plot routine output. 

The comparisons between the computer results and the exact solutions and 
test results are considered satisfactory. 
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CONCLUDING REMARKS 


An analysis and a computer program for predicting the transient thermal 
response of a charring ablation thermal protection system has been described. 
The numerical formulation of the equations is such that an implicit solution 
is obtained. This method of solution affords both a rapid and accurate solu- 
tion for both ablating and nonablating type problems . 

Provision is made in the program for a number of surface boundary condi- 
tions. These provisions allow efficient use of the program for analyzing both 
ground and flight test data and trajectory synthesis. 

The computer program has been checked out with both exact solutions and 
actual ablation test data. The numerical results are in good agreement with 
the exact solutions and test data. However, the analysis depends upon using 
good property values, and some effort must be expended in obtaining the best 
possible thermal properties. 


Manned Spacecraft Center 

National Aeronautics and Space Administration 
Houston, Texas, November 1 , 1965 


28 



REFERENCES 


1. Israel, Martin H. ; and Nardo, S. V.: An Annotated Bibliography on 

Ablation and Related Topics. PIBAL rep. no. 686 , Polytechnic 
Institute of Brooklyn, May 1964. 

2. Hurwicz, Henryk: Aerothermochemistry Studies in Ablation. Fifth AGARD 

Colloquium on High-Temperature Phenomena (Braunschweig, Germany, 

April 9-13? 1962, Pergamon Press, New York, 1963 . 

3. Liebmann, G. : A New Electrical Analogue Method for the Solution of 

Transient Heat Conduction Problems. Transactions of the ASME, vol. 

78 , no. 3 , 1956 , p. 655. 

4. Forsythe, George E.; and Wasow, Wolfgang R. : Finite-Difference Methods 

for Partial Differential Equations. John Wiley and Sons, Inc., 

New York, i 960 . 

5. Angelone, J., et al : Development of a Reinforced Carbonaceous and 

Ablative Composite for Entry Heat Protection of Manned Spacecraft. 

Tech. Rep. 330*17? Astronautics Div., Chance-Vought Corp., July 1963 . 

6 . Roberts, Leonard: Mass Transfer Cooling Near the Stagnation Point. 

NASA TR R -8 (Supersedes NACA TN 4391)? 1959* 

7. Dow, Marvin M. ; and Swann, Robert T.: Determination of Effects of 

Oxidation on Performance of Charring Ablators. NASA TR R-I 96 , 1964. 

8 . Carslaw, H. S.; and Jaeger, J. C.: Conduction of Heat in Solids. 

Clarendon Press, Oxford, 1959- 

9- Schneider, Paul J.: Conduction Heat Transfer. Addison-Wesley Publishing 

Co., Inc., c. 1955. 




APPENDIX A 


SAMPLE PROBLEM 


The following sample problem is shown to .indicate the form of the data 
input and the program output. A typical charring material subjected to a 
constant heating rate as experienced in an arc tunnel is presented. The fol- 
lowing is a sketch of the model: 

q = 95 Btu/ ft 2 - sec 
cw 

1 " 

1.5 in. 

_L 

0.1 in . 

Adi abaf ic 





The various material properties and dimensions are shown in the program 
output of Table- II. The insulation is assumed to be ablation material for this 
problem. The problem coding sheet and subsequent data card listing are shown 
in Table I. The initial temperature of the structure was assumed uniform and 
equal to 530 ° R ( 70 ° F). Figures 8, 9> and 10 are the output data obtained 
from the plot routines. 


31 



APPENDIX B 


PROGRAM USAGE INSTRUCTIONS 


IBM 709 4/40 program F021, standard ablation program, designated STAB II, 
is designed to evaluate the transient thermal performance of a charring ablation 
heat protection system. The program considers ojae ablating material and up to 
12 different materials in the supporting backup structure. A maximum of 
50 nodes may be considered in the ablation material, and a maximum of 10 nodes 
per material is allowed for each backup structure material. Air gaps can be 
considered between successive materials in the backup, thus allowing for 
both radiative and/or convective heat transfer between materials. The heat loss 
to the cabin environment from the backup structure can be accomplished by both 
radiation and/or convection, or an adiabatic backface surface may be prescribed. 

Unless otherwise specified, the input problem data are in "floating 
point" form (E12.8 format) and must end in columns 12, 24, 56, 48, 60, and 72. 

It is suggested that each floating point number have a sign, a two-digit 
exponent, and a decimal point. For example, the number 145- 25 can be written 
as +1.4525 +02, +145.25 +00, or +.14525 +05. 


Input Nomenclature 

The nomenclature used in the problem data input is as follows: 

NCASE number of problems to be run successively 

HEAD any 72 alphabetical and/or numerical characters 

TITLE control card for reading in new input for successive 
problems 

1. blank card — new data will be read in 

2. Six asterisks in columns 1 to 6. Skip to 

next read statement 

TLIM time limit of problem, sec 

TINT starting time of problem, sec 

NFTT number of points in time— step table (the minimum 

value of NPTT is 2) 

NPL/T output plot control 

=1 plot routine will be used 

=0 plot routine will be ignored 
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TTABLE 


DELTT 

IPRC 

FC0NV 

FRAD 

TABL 

TCHAR 

TREC 

RH0V 

RHjZ$C 

FBL0W 

EMV 

EMC 

H300 

VL 

HV 

VPT 


=0 

=1 

FV 

TV 


time in time— step table, sec 

time step to be used for each calculation — starting 
at time TTABLE, sec 

variable print frequency in TTABLE table; that is, if 
DELTT = 1. 0 and IPRC = 10, the output will be printed 
at 10— second intervals 

factor to correct convective heating rate for various 
body locations 

factor to correct radiative heating for various body 
locations 

temperature at which ablation starts, °R 

temperature at which ablation stops, °R 

surface temperature, °R, or time at which char removal 
is to start, sec 

density of virgin ablation material, lb/ft ^ 
density of mature char material, lb/ft ^ 

blowing efficiency of ablation gases in reducing convective 
heating 

emissivity of virgin ablation material 

emissivity of charred ablation material 

enthalpy of air at 300° K, 129* 0 6 Btu/lbm 

initial thickness of virgin ablation material, in. 

heat of degradation of virgin material, Btu/lbm 

test to determine if the reaction zone and char zone 
thermal properties are irreversible with temperature 

properties are irreversible and equal to the value at 
the maximum individual node temperature (this is the 
recommended value for VPT) 

properties are reversible 

view factor for external environment 

sink temperature of external environment, °R 
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CHARK -thermal conductivity of material at TCHAR, Btu/ft— hiy- 0 R 

CHARC specific heat of material at TCHAR, Btu/lbmr-°R 

ABLK thermal conductivity of material at TABL, Btu/ft— hr— °R 

ABLC specific heat of material at TABL, Btu/ lbirt-°R 

RP number of node points in ablation material 

RKC number of points in char thermal conductivity — temper- 

ature table 

RCPC number of points in char specific heat — temperature 

table 

RKV number of points in virgin thermal conductivity — 

temperature table 

BCPV number of points in virgin specific heat — temperature 

table 

RREC number of points in surface recession - temperature or 

time table 

TKC temperature values in char thermal conductivity — 

temperature table, °R 

XKC thermal conductivity values in char thermal conductivity 

temperature table, Btu/ft— hr— °R 

TCPC temperature values in char specific heat — temperature 

table, °R 

CPC specific heat values in char specific heat — temperature 

table, Btu/lbnt-°R 

TKV temperature values in virgin thermal conductivity — 

temperature table, °R 

XKV thermal conductivity values in virgin thermal conduc- 

tivity — temperature table, Btu/ft— hr— °R 

TCPV temperature values in virgin specific heat — temperature 

table, °R 

CPV specific heat values in virgin specific heat temperature 

table, Btu/lbm— °R 

TS temperature, °R, or time, sec, values in the surface 

recession table 



SR 

KTRAPT 

TIME 

QC01J 

QRAD 

VEL 

RMB 

NPBS 

BL 

XNPM 


surface recession values in the surface recession — 
temperature or time table, in. / sec 

number of time points in the trajectory input table 

the array of (NTRAPT) trajectory time values, sec 

the corresponding array of cold wall convective heating 
2 

rates, Btu/ft -sec 

the corresponding array of radiative heating rates, 

O 

Btu/ft - sec 

the corresponding array of flight velocity, ft/ sec 

number of materials in backup structure 

total number of node points in backup structure 

total thickness of backup structure, in. 

number of nodes in each individual material in backup 
structure 


KKPB number of points in each individual backup structure 

material thermal conductivity - temperature table 


RCPB number of points in each individual backup structure 

material specific heat — temperature table 

X3DNT any 72 alphanumeric characters used to describe each 
individual material in the backup structure 

TXK temperature values in backup material thermal conductivity 

temperature table, °R 

XK thermal conductivity values in backup material thermal 

conductivity — temperature table, Btu/ft— hr— °R 

TCP temperature values in backup material specific heat — 

temperature tables, °R 

CPX specific heat values in backup material specific heat — 

temperature tables, Btu/lbm-°R 

RH/BX density of individual materials in backup, lb/ft ^ 

XBM thickness of individual materials in backup, in. 

EMFB emissivity of front surface of each material in backup 
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EMBB 

H 

GAPX 

FTEST, 

BTEST 

=0 

=+l 

=-l 

TENV 

HEJW 

FENV 

Qh/SS 

=0 

=+l 

TEST2 

=0 

=-l 

=+l 

TEMPI 

TX0 


emissivity of back surface of each material in backup 

film coefficient between adjacent materials in backup, 
Btu/ft^-hr-°R 

width of gap between adjacent materials in backup, in. 

tests to determine the mode of heat transfer between 
materials for the front and backface of each material 
respectively 

conduction only between materials 
convective heat transfer only 

radiation only or radiation and convection heat transfer 

temperature of interior cabin environment, °R 

film coefficient to interior cabin environment, 
Btu/ft 2 -hr-°R 

view factor and emissivity product for radiative heat 
transfer to cabin interior 

boundary condition between last node of the backup 
structure and cabin environment 

adiabatic surfaces 

radiation and/or convective loss 

determines the proper heat shield initial temperature 
distribution 

constant, uniform initial temperature distribution 
arbitrary initial temperature distribution 
linear temperature distribution 

temperature to be used when constant temperature 
distribution option is used, R 

initial temperature at front surface of heat shield to 
be used in computing initial linear temperature 
gradient, R 

arbitrary temperature distribution values, to be used only 
if TEST2 is negative, °R 


TEMDI 



NHP number of points in enthalpy — temperature curve fit 

HX enthalpy values in enthalpy — temperature table, 

Btu/lbm 

TW corresponding temperature values in enthalpy — temperature 

table, R 


Input Data Card Preparation 

The input data are given in the following order. Each number in the fol- 
lowing listing refers to a separate record and must begin on a new data card. 
The input data have been grouped, where possible, into various sections dealing 
with a particular part of the input, that is, ablation material properties, 
trajectory data, backup structure, et cetera. This grouping permits the use of 
a minimum number of input cards for running successive problems. The title 
card as described in the input nomenclature controls the input for successive 
problems . 


1 . The first data card contains the value of NCASE. NCASE is an integer 
(15 format) and must end in column 5 . This card tells how many problems are to 
be run and is entered only once at the start of the data deck. 

2. Columns 1 to 72 of the second data card contain any title or identi- 
fication information desired; any alphanumeric character may be used. This 
card is printed at the top of the first page of the output. This card must be 
included in all successive problems to be run. 

(a) Problem time section 

3. TITLE card — if blank, cards 4 and 5 must be submitted. If six aster- 
isks are punched in columns 1 to 6 , skip to record number 6 . 

4. This record contains, in the following order, TLIM, TINT, NPTT, and 
NPL0T. TLIM and TINT are entered as floating-point numbers and must end in 
columns 12 and 24. NPTT and NPL0T are integers entered with an 15 format and 
must end in column 30 and 35 * 

5. Start entering the values of TTABLE , DELTT , IPRC. TTABLE and DELTT 
are floating-point numbers and must end in columns 12 and 24. IPRC is entered 
as integer with an 15 format and must end in column 30. Use as many cards as 
required to enter NPTT values. 

(b) Heating rate factors section 

6 . TITLE card - if blank, card 7 must be submitted. If six asterisks are 
punched in columns 1 to 6 , skip to record number 8 . 
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7. Enter the FC0NV and FRAD. These numbers are entered as floating-point 
numbers and must end in columns 12 and 24. 

(c) Ablation material section 

8 . TITLE card — If blank, cards 9 to l 8 must be submitted. If six aster- 
isks are punched in columns 1 to 6 , skip to record number 19 • 

9. HEADNG card — any alphanumeric characters in columns 1 to 72. Records 
9 to l 8 contain input data for the ablation material. 

10. Enter TABL, TCHAR, TREC, RH0V, RH0C, and FBL0W . These numbers are 
entered as floating-point numbers (6E12.8 format) and must end in columns 12, 

24, 36 , 48, 60 , and 72 . 

11. Enter EMV, EMC, H300, VL, HV, and VPT . Use the same format as card 10. 

12. Enter FV, TV, CHARK, CHARC, ABLK, and ABLC . Use the same format as 
card 10 . 


13. This card contains, in the following order, UP, NKC, NCPC, NKV, NCPV, 
and NREC. These numbers are fixed-point integers and must end in columns 5> 

10, 15, 20, 25 , and 30. An 15 format is used to read in these numbers. 

14. Start entering the curve of TKC versus XKC, with the values of TKC 
ending in columns 12, 36 , and 60. The corresponding values of XKC must end in 
columns 24, 48, and 72; for example, three TKC-XKC points are contained on one 
card. The numbers are entered as floating-point numbers. Use as many cards as 
required to enter NKC points on the curve. 

15 . Start entering the curve of TCPC versus CPC with the values of TCPC, 
ending in columns 12, 36 , and 60. The corresponding values of CPC must end 

in columns 24, 48, and 72; for example, three TCPC-CPC points are contained on 
one card. The numbers are entered as floating-point numbers . Use as many 
cards as required to enter the NCPC points on the curve. 

16. Start entering the curve of TKV versus XKV with the values of TKV 
ending in columns 12, 36 , and 60. The corresponding values of XKV must end in 
columns 24, 48, and 72; for example, three TKV -XKV points are contained on one 
card. The numbers are entered as floating point. Use as many cards as re- 
quired to enter the NKV points on the curve. 

17. Start entering the curve of TCPV versus CPV with the values of TCPV, 
ending in columns 12, 3 6 , and 60. The corresponding values of CPV must end in 
columns 24, 48, and 72; for example, three TCPV-CPV points are contained on one 
card. The numbers are entered as floating point. Use as many cards as re- 
quired to enter NCPV points on the curve. 

18. Start entering the curve of TS versus SR with the values of TX, ending 
in columns 12, 36 , and 60. The corresponding values of SR must end in columns 
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24, 48, and 72; for example, three TS-SR points are contained on one card. 

The numbers are entered as floating point. Use as many cards as required to 
enter NREC points on the curve. 

(d) Trajectory data section 

19. TITLE card - if blank, cards 20 to 22 must be submitted; if six 
asterisks are punched in columns 1 to 6, skip to record number 23. 

20. HEADNG card - any alphanumeric characters in columns 1 to 72. Records 
21 and 22 contain trajectory input data. 

21. Enter ETRAPT . This number is an integer and must end in column 5* 

An 15 format is used to read in this number. 

22. Start entering the trajectory data in the following order: TIME, 

QC0E, QRAL? VEL . These values are entered as floating-point numbers and must 
end in columns 12, 24, 36, and 48. There are four trajectory data points on 
one card. Use as many cards as required to enter ETRAPT points in the tra- 
jectory. 


(e) Backup structure section 

23. TITLE card - if blank, cards 24 to 31 must be submitted; if six 
asterisks are punched in columns 1 to 6, skip to record number 32. 

24. Enter EMB, EPBS, and BL. These three values must end in columns 5? 
10, and 24. EMB and EPBS are integers and are read in under an 15 format. BL 
is a floating-point number. 

25* Enter the values of XEPM. XEPM is in floating-point form and must 
end in columns 12, 24, 36, 48, 60, and 72. Use as many cards as required to 
enter EMB points. 

26. Enter the values of EKPB and ECPB 0 These numbers are integers and 
EKPB must end in columns 5? 15? 25, 35? and 45; and the corresponding values of 
NCPB must end in columns 10, 20, 30, 40, and 50. An 15 format is used to read 
these values. Five EKPB-ECPB values are contained on one card. Use as many 
cards as are required to enter EMB points. 

27. XIDET card - any alphanumeric characters in columns 1 to 72 . This 
card contains a description of each backup material. 

28. Start entering the curve of TXK versus XK with the values of TXK, 
ending in columns 12, 36, and 60 o The corresponding values of XK must end in 
columns 24, 48, and 72; for example, three TXK-XK points are contained on one 
card. The numbers are entered as floating point. Use as many cards as re- 
quired to enter EKPB points on the curve. 
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29. Start entering the curve of TCP versus CPX with the values of TCP, 
ending in columns 12, 36 , and 60. The corresponding values of CPX must end 
in columns 24, 48, and 72; for example, three TCP-CPX points are contained 
on one card. The numbers are entered as floating point. Use as many cards as 
required to enter RCPB points on the curve. Repeat records 27, 28, and 29 
until the properties for 3MB materials have "been entered. The maximum number 
for 3MB is 12. 

30. Start entering the following values in order: RH0BX, XBM, EMFB, and 

EMBB. These values are entered as floating-point numbers (6E12.8 format) 
and must end in columns 12, 24, 36 , 48, 60, and 72. Use as many cards as 
required to enter 3MB points. 

31. Start entering the following values in order: H, GAPX, FTEST, and 

BTEST. These values are entered as floating-point numbers (6E12.8 format) 
and must end in columns 12, 24, 36 , 48, 60 , and 72. Use as many cards as 
required to enter 3MB points. 

(f) Interior environment section 

32. TITLE card - if blank, cards 33 and 34 must be submitted; if six 
asterisks are punched in columns 1 to 6 , skip to record number 35 * 

33* HEADUG card — any alphanumeric characters in columns 1 to 72. 

Record 35 contains properties of environment. 

34. Enter the following: TE3W, HE3)TV, FENV, and QL0SS. The values are 

entered as floating-point numbers and must end in columns 12, 24, 36 , and 48. 

(g) Initial temperature section 

35. TITIE card - if blank, records 36 and 37 must be submitted; if six 
asterisks are punched in columns 1 to 6, skip to record 39. 

36 . HEADNG card — any alphanumeric characters in columns 1 to 72. 

Records 37 and 38 contain initial temperature distribution input. 

37- Enter TEST2, TEMPI, and TX0. These values are entered as floating- 
point numbers and must end in columns 12, 24, and 36 . 

ROTE: If TEST2 is a negative number, record 38 must be submitted; other- 

wise, skip to record 39 ■- 

38 . Enter the arbitrary temperature distribution values, TEMDI. These 
values are entered as floating points with a 6E12 . 8 format. Use as many 
cards as required to enter NP plus NFBS node points. 

(h) Enthalpy - temperature section 
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39- TITUS card — if blank, records 40 and 4l must be submitted; if six 
asterisks are punched in columns 1 to 6, this is the last data card in the 
problem input. 

40. Enter NHP. This value is an integer and must end in column 5* An 
15 format is used to read in this number. 

41. Start entering the curve of HX versus TW with the value of HX ending 
in columns 12, 56 , and 60. The corresponding values of TW must end in columns 
24, 48, and J2; for example, three KX-TW points are contained on one card. 

The numbers are entered as floating points. Use as many cards as required to 
enter KEP points on the curve. Record 4l consists of the last data cards 
required as input for a problem. 

As many successive problems as desired may be run at one time by proper 
input preparation. STAB II has been designed to save all input information 
until it is changed by new input data. Therefore, the use of the TITLE 
control card is very important when running more than one problem and using 
the input data of previous problems. As shown, each input section starts 
with a TITLE control card for determining whether new input data are to be 
used. If any data are changed within a section, then all data cards required 
for that section must be submitted. 

STAB II can also be used for solving one— dimensional transient heat- 
conduction problems of nonablating materials. The following input parameters 
must be adhered to: 

1. TABL must be greater than the maximum temperature expected during the 
calculation. Also, TABL )> TCEAR )> TREC. 

2. The ablation material must be considered to be the first material 
in the structure for calculation purposes. 

3 . The virgin and char properties must be inputed as described above 
but can have the same values; that is, XKV = XKC, CPC = CPV, RH0V = RH0C, 
et cetera. 

The following dimensional statements and program limitations should not 
be violated when preparing the input described above for ablating and non- 
ablating structure: 

1. All property tables can have a maximum of 20 points (i. e. , a tempera^ 
ture and specific heat value constitute one point). 

2. The surface recession table can have a maximum of 50 points (TS and 
SR constitute one point). 

5- The trajectory table can have a maximum of 300 points (TIME, QC0N, 
QRAD, and VEL constitute one point). 
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4. The ablation material can he broken into a maximum of 50 nodes. The 
backup structure can consist of up to 12 different materials with a maximum of 
10 nodes per material. 

5. A minimum of three nodes per material (ablation or backup) must be 
specified. 

6. A minimuim of two materials must be specified (ablation material and 
one backup structure material). 

7. Pure conduction only is allowed between the ablation material and the 
first material in the backup. 

8. If any data input is changed in the Ablation Material Section on 
successive problems, the Ablation Material Section data cards plus the 
Initial Temperature Section data cards must be submitted. 


Program Output Information 

The computed results are available in two forms of output: tabular and 

plot outputs. The tabular output presents the computed results in block 
type form for each computation step as controlled by the print count control 
number. As discussed in the preparation of input data, both the computational 
time step and print control can be varied throughout the running of a problem. 
Therefore, excessive printed output is avoided, and there is a considerable 
savings in actual machine computation time. The plot outputs are printed 
and plotted only when the entire set of problems to be run are completed. 

Tabular output . — The program prints a listing of the data input para- 
meters for identification of the problem and ease in identifying any input 
mistakes. For stacked problems, the program prints only that input information 
that is changed from the previous problem. The following calculated problem 
output is printed: 

1. Time, sec 

o 

2. Cold wall convective heating rate without blowing, Rtu/ft -sec 

3. Radiative heating rate, Btu/ft^-sec 

4. Velocity, ft/sec 

5. Gas a lation rate, lbm/ft^-hr 

6. Char ablation rate, lbm/ft -hr 

7. Total ablation rate, lbm/ft 1 2 3 4 5 6 7 8 -hr 

8. Surface recession depth from original surface, in. 


^3 


2 

9* Hot wall convective heating rate without blowing; Btu/ft -sec 

10. Temperature distribution in ablation material; °H 

11. Temperature distribution in backup structure; °R 

The temperatures printed for the ablation material are for fixed distances 
from the original surface. These distances are calculated from the initial 
ablation material thickness and number of nodes in ablation material. For ex- 
ample 

let 


VL = 1.0 in. 


HP = Id- 


then 


AX = 



0.1 


The temperatures will be printed for X distances of 0; 0.1; 0.2; 0-3* 
et cetera; from the original surface until the surface has receded beyond these 
fixed distances at which time the node no longer exists and is dropped from 
the printout. This is illustrated in the following way: let surface 

recession = 0.2o inch. The first temperature printed then is the surface 
temperature of the material; located 0. 26 inch from the original material sur- 
face. The following printed ablation material temperatures are for X dis- 
tances of 0. 0.4; 0.5; ...; 1.0 inch. 

The format for the temperature distribution printout is El6. 5 with six 
temperatures printed per line. 

Plot output . — The plot output gives the following ablative material per- 
formance parameters as a function of time: 

1. Surface depth; in. 

2. Bondline temperature between ablator and backup structure; °R 

3. Two selected isotherm depths 

These values are also printed in tabular form for ease in checking and 
replotting of the results. The plotted curves contain all maximum and mini- 
mum values of the parameters. 
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IBFTC MAIN AOOnn 

Aooin 

STRUCTURES AND MECHANICS DIVISION A0020 

THERMO-STRUCTURES BRANCH A0030 

THERMAL PROTECTION SYSTEMS SECTION AOOun 

A005fl 

THIS PROGRAM DETERMINES THE PERFORMANCE OF A CHARRING ABLATOR A0060 

A0070 

ANALYSIS AND PROGRAM DEVELOPED By DONALD M. CURRY * ES32 A0080 

A0090 

DIMENSION ESAVE1 (3) »ESAVF?<3) ,ESAVE3(3) A0100 

DIMENSION TITLE ( 12 ) , HFADNG ( 12 ) , XlDNT (12.12) »TKC ( 20 ) » XKC < 20 ) » A0110 

1CPC (20 ) , TK V ( 20 ) *XKV(20) ,TCPV(20 >, CPV ( 20 >. TIME < 30(1 ) .OCONOOO) , AO 120 

20RAD ( 300 ) , VEL ( 300 ) , XnPM ( I 2 ) .NKPB ( 12 ) , NCPR ( 12 ) ,TXk ( 20 . 12 ) , XK ( 20 , 12 ) A0130 

3, TCP ( 20 , 12) , CPX (20,12) * RHOBX ( 12 ) , XBM( J 2 ) .FMFB(1?) »EMBR<12) ,HXX< 12) AOlttO 

4 , GAPX ( 12 ) , FTFST ( 12) »BTEST (12) .TFmDI ( 200 ) » TX1 (200) .TX2(200) , AO 150 

5TXPT< 10* 12) ,TUL1 (200) , TUL2(200) » HX (50 ) , Tw ( 50 ) , IP ( 50 ) , IR1 ( 50 ) , AO 160 

6lP2<50) »TUL<5(1> » IFM( 50 ) . TY ( 200 ) , A ( 200 ) ,B(200> »C(200) ,D<200> , A0170 

7R ( 50 ) * RHO (50 ) ,CP(50) ,DXB(12) ,XKR(10,12) ,CPB(10,12> .XMDG(50) , A0180 

flYK (50) ,AB(10,l2),BB(in,12).CR(10,12),DB(10»12)»SR(10,12), AO 190 

9RR1I10.12) ,RB2< 10,12) , H < 1 2 ) »S (50 ) , NPM ( 12 ) A0200 

DIMENSION TTUL ( 50 ) * RHOY1 ( 50 ) ,PHOy2 < 50 ) .DRHO ( 50 ) * TCPC ( 20 ) A0210 

DIMENSION TlMFP ( 300 ) , PRFs ( 300 ) * Xc ( 50 ) ,TX2C(50) .XV (50) ,XDV(50) A0220 

DIMENSION TS ( 50 ) . SR ( 50 ) A0230 

DIMENSION TT ARLE ( 20 ) . DELTT (20) .IPRC(20) A0240 

DIMENSION ASAVE1 (3) . ASAVF2 ( 3 ) »ASaVE3(3) .BSAVE1 (3) .BSAVE2(3) . A0250 

1RSAVE3(3> »CSAVE1(3) .CSAVE2I3) .CSAVE3(3) .HFAD( 12) . A0260 

1DSAVE1 (3) »DSAVE2(3) »DSAVF3(3) A0270 

DIMENSION XRA(30) ,YA(30) A0280 

C A029O 

COMMON TKC. XkC.TcPC. CPC. TKV.XKV.TCPV.CPV.XNPM. RHOBX. XRM.FMBR, A0300 

lFMFB.NKPb.NCPP.TXK.XK. TCP . CPX ,NPM. GAPX .FTFST . RTEsT « TEMDI , TX1 , A0310 

2TX2.TX2T.TUL.TULl.TUL2. IR.IR1.IR2. A.B.C.D.S.R. AP.BB.CR.DB.SB, A03PO 

3RR1.RB2.TY. RHOY1 , RH0Y2. XMDG.RHO, CP » YK .XKR.CPB.DXp.DT, XLOST . A0330 

4T ABL .TCHAR.TRFC .RHOV . RHOC .FBLOW . fMV . EMC .H300 , NKC, NCPC .NK V. NCPV . A0340 

5NP . NMR . NPRS . NPF .TEST 2 . TEMPI .TXO » TFNV. HENv . FFNV. OLOSS , TLIm , TINT A0350 

COMMON II. 12. 13. 14. 15. 16. OIN. I NT. DX. XMT .TL » VL . BL , DMP, ERR 1 . ERR2 » A0360 

1fPR3.FRR4.HV, VPT.CHAPK.CHARC. ABLk. ABLC.XmDC.H A0370 

C A0380 

3000 FORMAT ( 12A6) A0390 

3001 F0RMAT(1X,12A6) A0400 

3002 FORMAT (6E12.8) A0410 

3003 FORMAT ( 615 ) A0420 

3004 FORMAT (115) A0430 

3005 FORMAT ( 215 ) A0440 

3007 FORMAT (2I5.lE14.fl) A0450 

3008 FORMAT (///1X.12A6) A0460 

3009 FORMAT (1H1.1X.12A6) A0470 

3010 FORMAT (4E12. 8 ) A0480 

3011 FORMAT (2E12.8.I6. I5.F13.8.E12.8) A0490 

3012 F0RMAT(2E12.8.I6> A0500 

DATA PRV0US/0545454545454/ A0510 

RFWIND 11 A0520 

ST0P=9999. A0530 

READ (5,3003)NCASE A0540 

LPLOTsO A0550 

JCNT=0 A0560 
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SO NK=1 A0570 

11=2 A05S0 

12=2 A0590 

T 3=2 A0600 

14=2 A0610 

15=2 A0620 

15=2 A0630 

117=2 A0640 

INT=1 A0650 

XLOST =0.0 A06S0 

XMT=0 . 0 A0670 

XMDT=0.0 A06R0 

FRR1=0.0 A0690 

FRR2=0.0 A0700 

FRR3S0.0 A0710 

FRR4=0.0 A0720 

ICT=0 A0730 

ICONT=0 A0740 

XMDCSO.0 A0750 

NKP=1 A0760 

XL STV=0 . 0 A0770 

NRS=2 A0780 

FRR5=0.0 A0790 

tpct=o AOftnn 

ICTP=0 A0810 

IPL0T=1 A0820 

NXA=1 A0830 

NXB=1 A0840 

NXC=1 A0850 

NXD=1 AOBSO 

BAVY3=-1()0. A0870 

SAVY4XS-100. A0880 

5X0=0. 0 A0890 

5DOT=0.0 A0900 

C A0910 

C GFNERAL TITLF OF PROBLEM A0920 

100 RE AD (5.3000) (HEAD ( K > , K=1 , 12 ) A0930 

WRITE (6 » 3009) ( HEAD ( K ) , K=1 . 12 ) A0940 

LPLOT=LPLOT+l A0950 

WRITE < 11 > <HEAD( I ) . 1=1. 12) A0960 

WRITE ( 6 . 1 1 0 ) A0970 

1)0 F0RMAT(//1X.11HINPUT DATA.//) A09B0 

READ(S.3000) (TITLE(L) »L=1 .12) A0990 

TF(TITLE(1J .F0.PRV0U5) 60 TO 150 A1000 

RFAD(S»3011)TLIM.TINt.NPTT,NPLOT,DMP,TDMp A1010 

RFADI5.3012) (TTABLE(I).DFLTT(I),IPRC(I),I=1.NPTt) A1020 

T=TINT A1030 

DTS=DFLTT(1) A1040 

DT=DELTT (1)/ 35 00.0 A10R0 

WRITE ( 6 . 120 ) TLIM.TINT.MPTT A1050 

120 FORMAT(lHO. 11HTIME LIMIT=» 1PF1 0.4.4X.13H INITIAL TIME=, 1PF1 0 . 4 . 4X , 5 A1070 

1HNPTT=»I4) A1080 

WRITE (6 . 122 ) A1090 

122 FORMAT ( //8X » 4HTIME . 1 OX . 9HTIME STfP . 6X , 13MPRI NT CONTROL) AllOO 

WRITE (6 . 124 > ( TT ABLE ( T ) . DFLTT (I) ,IPRC(I) .I = 1»NPTT) A1110 

124 FORMAT (5X ► 1PE1 0.4.6X, ) PF10 .4 , 9X . 14 ) A1120 

C A1130 
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c location factors for CONVFCTIVE and RAOIATIVE heating 
150 READ(5* 3000 ) ( TITLE (l) * L=1 * 12 > 

IF(TITLE(1> .EO.PRVOUS) GO TO 200 
READ(5*3002) FCONVtFRAO 
WRITE (6* 155 ) FCONV*FRAD 

155 FORMAT ( 1H0 * 6HFC0NV= » 1PE12 . 5 * 4X5HpR AD= * 1PF 1 2 . 5/ ) 

PROPERTIES OF ABLATION MATERIAL 
200 RFAD(5,3000) (TITLE(L) *L=1 * 12 ) 

IF (TITLE 111 .EO.PRV0U5) GO TO 300 
READ ( 5» 3000 ) (HEADNG(K) ,Krl.l2) 

RF AD ( 5 * 3002 ) T ABL * TCH AR. TREC , PHOv , RHOC * FRLOW , FMV , EMC * H300 , VL , HV * 

1 VPT * F V * T V * CHAPK , CHAR C , ABLK * ABL C 
RFA (5*3003) NP*NKC*NCPC,WKV,NCPV*NREC 
READ (5* 3002) (TKC(K)t XKC (K ) *K=1 * NKC ) 

READ(5*3002) (TCPC(M) , CPC (M > ,M=1 , NCPC ) 

RFAD(5*3002) { TKV (L ) , XKV <1. ) *L=1 *NKV ) 

RFAD (5*3002 ) (TCPV(N) ,CPV(N) ,N=1,NCPV) 

Rf AD(5*3002) (TS(I) ,9R(I > *I=1,NRFC) 

WRITE(6*3008) (HEADNg ( K ) ,K=1 . 12 ) 

WRITE (6*210) TABL*TChAR* TREC , RHO v * RHOC . FRLOW , FM V , EMC , H30 0 * VL , HV » 

1 VPT , F V * T V . CHAPK . CHARC * ABLK , ABLC 

210 FORMAT ( 1H0 *5HTABL= *1PF12.5*3X ,6HTCHARs» 1PE12.5*3x* SHTREC s* 1PE1 2,5* 
13X »5HRH0v = * 1PF12 .5* 3x* 5HRH0Cr * 1PF12.5* 21 X/1X, 6HFBL0W=, 1 PE 12. 5* 4X, 4 
2HEMV=,1PE12.5*4X,4 HEmC=*1PE12.5»3X*5HH300=*1PE12*5.5X,3HvL=. 1PE1?, 
35/i»X,3HHV=*lPFl2.5*4x*4HVPT=.lPEl2.5,5X»3HFV=*lPFl2.5,5X,3HTVr,lPE 
1 12. 5* 2X*6HCHAPK=,lPEl2. 5/1 X*6HCHARC=*1PE1 2. 5* 3X*5HABLK=*1P£12.5*3X 
2,5HABLC=*lPEl2.5/) 

VLISVL 

VL=VL/12.0 

VLV=VL 

WRITE (6*220) nP.NKC*NCPC.NKV.NCPv*NREC 

220 F0RMAT(2X.3HNP=*1I4,4X*4HNKC=,1I4*4X,5HNCPC=,1I4,4X*4HNKV=.1I4*4X. 
1 BHNCPVS * 1 14 » 4X * 5HNRECS *114) 

WRITE (6*221 ) 

221 FORMAT (/32XH5HVIRGIN MATFRIAL/2nX»7HTHERMAL.3BX,8HSPFClFlC/3X. 11 H 
1TFMPERATURE*4X,12HC0N0UCTTVITY.19X.11HTEMPERATURF*7X,4HHEAT) 

KLLL=MIN0 ( NK V . NCPV ) 

WRITE (6*222) ( TK V (L) , XKV (L ) *TCPV (L ) *CPV (L) *L=1 *KlLL) 

222 FORMAT (2X * 1PE! 2.5* 4X, 1PE12.5*18X, 1PE12 *5 *3X,lPEi2.5) 

IF(NKV-NCPV) 223,227,225 

223 kLLLL=KLLL+1 

WRITE (6 * 224 ) (TCPV(L) ,CPV (L > .LrKLLLL, NCPV) 

224 FORMAT (48X * 1PF12.5* 3x , 1PF1 2.5) 

GO TO 227 

225 kLLLL=KLLL+1 

WRITE (6 *226) ( TKV ( L ) , XKV (L ) *L=KLlLL *NK V) 

226 FORMAT ( 2X , 1PE1 2,5*4X,lPEt2,5) 

227 WRITE(6*228) 

228 FORMAT (//33X,14HCHAR MATFRI AL/20X *7HTHERMAL. 38X, 8HSPFCIFIC/3X, 11 H 
1TFMPERATURE * 4X * 12HC0NDUCT I V IT Y , 1 9X * 1 IHTEmPER ATURf * 7X , 4HHE AT ) 

KLLLrMINO ( NKC » NCPC ) 

WRITE (6* 222) ( TKC (L ) , XKC (L ) *TCPC(L) *CPC (L ) *L=1 *KlLL) 

IF(NKC-NCPC) 230,235,232 
230 KLLLL=KLLL+1 

WRITE (6* 224 ) (TCPC(L) * CPC (L) *L=KLLLL* NCPC ) 

GO TO 235 
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on on 


232 Kt LLL=KLL! +1 A1710 

wP I TE ( 6 » 226 ) <TKC(L) ,XKC(I. > »L=KL| LL.NKC) A1720 

235 WRITE ( 6 » 240 ) A1730 

240 FORMAT(//28X.23HSIIRFaCE RFCFSFION T ABLE/ /25X , 1 IHtFMPFR ATi |RE , AX , 1 1 H A1740 

ISP - TN/SFC) A17F0 

WR I TE ( 6 *245 ) (T5(T) rSP(T) rI=l.NPFC) A17A0 

245 FORMAT ( 24X . 1PF12 . 5 , 7x , 1PF 12.5) A1770 

A1780 

PROPERTIES OF TRAJECTORY A1790 

300 RFAD(F,3u00) (TITLEIL) »L = 1 .1?) AIAOO 

lF(TITLE(1 ) .FO.PRVONc;) PO TO 400 A1R10 

RFAD(F,3000) (HFAONPIL) >1 =lfl2> A1A?0 

RFAO(F.3U04) MTPAPT A1R30 

RFAD(F»3010) (TTmF(K) .QCONIK ) ,QRAD(K) , VFL (K) .K = l .NTRAPT) A1A40 

WRITE (6» 3008) ( HE ADN<; ( L ) , I =1 • 1 2 ) AlflSO 

WRITE(6»3I0) MTPAPT A1R60 

310 FORMAT! 1H0»27H NO. OF TRAJECTORY POINTS =rll4) A1R70 

wPI TE ( 6 » 320 ) A1A80 

320 FORMAT (//PX.APTIMF > A X » 1 ?HO CONVFcT I VE 1 4X , 1 1HO R An! ATI VE > 7X , RHVFl OC A1R90 
1 TTY ) A 191)0 

WRITE I 6» 330) ( T IMF (K ) , QCON(K ) ,ORaD(K) , VFI.(K) ,4 = 1 .NTRAPT) A 1910 

330 F0RMAT(lP4Elft.5) A1920 

A1930 

PROPERTIES OF BACK-Ilp STRUCTUPE A1940 

400 KF AD ( F » 3(100 ) (TITLE(L) *L=1 .12) A1QF0 

IF (TITLF! 1 ) .FO.PRVOIK) Ro TO F00 A19A0 

WRITE(6»410) A1970 

410 FORMAT(//10X,3lH PRORFRTTFS OF RACKUP STRUCTURE/) A1980 

RF AD ( F » 3ll 07 > mMP.NPRr.BL A1990 

kFAniF,3U02) (XNPM(K) ,K=1 .NMR) A2000 

RF AD ( F , 4 1 F ) (MKPH(I) .NCPR(I) .1=1, NMR) A2010 

415 F OKM AT ( 1 u T 5 ) A2020 

DO 420 K=1,Nmp A2030 

NPM(K)=XNPM(K)+0.000OO0O2 A2040 

420 CONTINUE A2050 

wPITE(6r42b) NMP.NPRq.Bt A20RO 

4?5 F0RMAT(/4X»3FHN0. OF MATFRIALF IN BACK-UP SHI FLDr • 1 I4/4X , 40HTOT At. A2070 

1NUMRER OF NODFS IN RaCK-IIP SHIELn=*lI4/4X*2RHTHICKNFSS OF BACK-UP A2080 
2sHIFLO=.lPEl2.5//) A2090 

RL=PL/1?.0 A2100 

no 440 1=1 »NMR A21 l 0 

LKrNKPBU) A2120 

lCP=NCPRlT) A2130 

RFAD l F t 3u00 ) ( ( X I nNT t K . I ) ) * K = ) , 1 ? ) A2140 

RF AD ( F » 3li02 ) ( (TXK ( J, I ) ,XK ( J, I > 1 , J=1 ,LK) A21F0 

RFAD(F.3ii02) ( (TCP(J.T) »CPX(J,I) ) , J=1,LCP) A21A0 

wP I TE ( 6 * 4 32 ) <XlDNT(Kt I) ,K = 1,12) A2170 

432 FORMAT ( / / 1 2A(S ) A21AO 

wP I TE ( 6 * 433 ) A2190 

433 F OR MAT ( //20X , 7HTHFRM aL » 3A Y t RHFPFCIF IC/3X . 1 1HTFMPFR ATtlRE . 4X , 1 ?HCOND A220O 

lllCTlVTTY. 19X, 1 1HTFMPFP ATtlRE .7X.4HHE AT) A2210 

KL LLzMINO (UK ,1 CP) A2220 

DO 434 N=),KLIL A2230 

WR 1 TE ( 6 . 222 ) (TXk (N»t ) ,XK (N. I ) ,TcP(N, I ) ,r.PX(N, I) ) A2240 

434 CONTINUE A2250 

IF(LK-LCP) 43F r 440 1 437 A22AO 

435 KLLLL=KLLl+l A2270 
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no 43 6 NzKLLLI »LCP 

WR I TE ( 6 r 224 ) (TCP (N» T > » CPY ( N r T ) ) 

436 CONTINUE 
GO TO 44 0 

437 kLLLL=KLLL+1 

no 43ft N=KLI_LL*Lk 

WRITE <6 r 226) (TXK (Nr T ) rXK (Nr I ) ) 

438 CONTINUE 
440 CONTINUE 

RF AD ( ft r 3n 02 ) (RHORX(L) r XBM(L) , EM fF ( L ) ,EMrB(L) rL=t r NMr ) 

RFAD <5# 3002) (H( J) ,OAPX ( J) r ETEsT ( J) r BTFST ( J ) r JrtrNMR) 

wP ITE ( 6 # 4S0 ) 

4B0 EORMAT(///55Xr 10HEMIsBlVTTY/8XrftHMATEPIALt 5X.7H0FNSlTYr7Xr9HTHICKN 
1FSS • 7V r ShFRONT r 9X r 4HRACK »7X r 14HNOPES/MATFP I Al_/ ) 

00 460 LLJ=1#NMR 

WRITE (6r4Sb) LLJr RHO rx ( L L J) rXPM(LLJ) r EMF r(LLJ) #FmRB(LLJ) » XNPM (LL j ) 
4S5 FORMAT ( UXr 111 r 8X r 1PF 1 0 . 4 #4X r 1 PEl 0 • 4 r 4X r 1 PE 1 0 . 4 r 4X r 1PF1 0 # 4 r 6X r 1PF1 
10*4/) 

460 CONTINUF 

wRITE(6f 465) 

465 FORMAT ( //4X r 60HADDI T I ONA( DATA FOR INDIVIDUAL MAtERIAI S TN BACKUP 
1 STRUCTURE // 1 l Y r ftHMATFR I At r 5X r 1 6HFILM COFFF IC I FNT , SX r 1 3HG AP THTCKnF 
2sS r 8X r 5HR TEST # 13X r 5HHTEST ) 

DO 480 Jr 1 9 NM° 

wRITE(6#470) Jr H ( J ) r GAPX ( J ) , FTFST ( J) r RTFST ( J ) 

470 FORMAT( l3Xr lT3 f l2Xr lPFln # Ur9XrlPFlO .4r 7X, 1PF1 1 . 4 , 7X r IPE1 1 . 4/ ) 

4ft0 CONTINUF 

PROPERTIES OF ENVIRONMENT 
500 RFAD ( Sr 3000 ) ( T I TLE ( L ) # L = 1 r 1 2 ) 

tF(TITLE( 1 ) .FO.PRVOUS) 00 TO 600 
RFAD ( Sr 3000 ) ( HEAONG ( L ) ,|„rlrl?) 

KF AD ( S r 3 u 02 ) TENVr HENVrFFNVr 01 OSS 
WRITE (6 #3008) (HFADN(;(L) ,1=1 #12) 

WRITE (6 r 520) TENV # HFi\lV# FFNV # QL OSs 

520 FORMAT ( /4X r 1 ?MTEMPERaTUPF=# 1PF12 # Sr 4X# 17HFILM COfFFIC TENT= # t PF12. S 
1 » 4X r 1 2HV I FW F A CTOR= # 1 PE 1 2 . 5 # 4X » 7HG L0ST=. 1PF12.5) 

INITIAL TFMPEPATUPE OTSTRTBUTTON 
600 RFAD ( Sr 3000 ) ( TITLE (L )# L = 1 # 12 ) 

TF(TITLE( 1 ) .FO.PRVOUB) Go TO 700 

RFAD (Sr 300 0) (HFAONG(L) rLrlrl?) 

j\jPF=NP+NPPS 

TLrVL+BL 

xNPrNP 

nX=VL/(XNP- 1 . 0 ) 

OXX=DX 

RF AD (Sr 3f)02 ) TEST? r TfmPI , TXO 
IF ( TEST2 ) 61 0 r 620 r 620 
610 RF AD ( S r 3u 02 ) ( TEMPI ( « ) r K = 1 r NPF ) 

HO 61 S KrlrNPF 
TX1 (K )=TEMDI (K ) 

TX2(K)=TX1 (K) 

TUL1 (K)rTXi(K) 

TUL2 (K)rTXl(K) 
b 1 5 CONTINUF 
l.=NP+1 
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no 6io 1 = 1 »nmp 

I N = NPM(I) 
no 617 J=1 r LN 
TX2T(JrI) =TFMni ( L ) 

I =L+1 

617 CONTINUE 
b 1 9 CONTINUF 
00 TO 626 
620 CALL TEMPO 

625 wR I TE ( 6 r 3008 ) ( HE ADN fi ( L ) , 1 =1 r 1 2 ) 

TF ( TEST2 ) 630»636#64o 
630 WRITE ( 6 9 632 ) 

632 FORMAT (4Xr 52HTEMPFRATI tRF niSTP I BuTlON IN HEAT SHlFLD IS APBRITARV/ 
1 ) 

wPlTE(6f 633) (TFMOI <K> r K = 1 rNPF) 

633 FORMAT ( 1P8E12 . 5 ) 

GO TO 64b 

635 wP I TE ( 6 r 637 ) TEMPI 

637 FORMAT ( //4X r 64HTEMPFR ATURF DISTR rFUT I ON TN HEAT SHIFLD Is UNIFORM 
lAND EOUAL TO » 1PE1 0 . 4/ ) 

GO TO 64b 

640 WRITE ( 6 9 64 1 ) 

641 p ORMAT ( 4X r 54 HI I NEAR TFMPfRATURE DISTRIBUTION ASSUMED IN HFAT SHIfL 
in/) 

WRITE (69 633) (TEMOI (L ) * L = 1 #NPF) 

645 IF(DMP) 7D0 r 700 * 646 

646 wR I TE ( 6 * 647 ) 

647 F0RMAT(//> 

648 wPITE<6#649> < TX 1 ( L ) . TX2 ( l ) r L=1 r NPF > 

649 F0RMAT(2XrlPFl2.br4X, 1 PF 1 2 . 5 ) 

WF I TE ( 6 r 650 ) 

650 F ORMAT ( // ) 


ENTHALPY AS A FUNCTION OF TFMPERATUPE 
700 «F AD ( s r 3n 0 0 ) (TITLE(L) »L=1 1 12) 

tF<TITLF ( 1) .EO.PRVOUS) GO TO 725 
RF AD ( 5 r 3004 ) NHP 

READ (Sr 3002) < HX < K > rTW(K) f K = l f NHp) 

725 DO 728 1=1 rNR 
" T P ( I ) =0 

T P 1 ( I )=0 
TP2 ( I )=0 
TFM< 1 ) =0 

xmdg( t )= n .o 

728 CONTINUF 

wPlTE(6r 730) 

730 FORMAT ( 1H1 r 12HOUTPUT DATA.//) 

XC( 1 )=0.0 
no 740 I =2 r Np 
XC(I)=XC(T-1)+DX 
740 CONTINUE 

7S0 TF ( T-TIME (NK ) ) 765r770r760 
760 nK=NK+1 

TF(NK-NTkAPT) 750 #750,762 

762 WRITE ( 6 r 763 ) MK 

763 FORMAT ( 1 HO r 3.3H THF VALUF OF NK 
GO TO 905 


A2850 
A2B60 
A2870 
A2880 
A2890 
A2900 
A2910 
A2920 
A2930 
A2940 
A2950 
A2960 
A2970 
A2980 
A2990 
A3000 
A3010 
A30P0 
A3030 
A3040 
A3050 
A 3060 
A3070 
A3080 
A3090 
A3100 
A31 10 
A31 ?0 
A 31 30 
A3140 
A3150 
A3160 
A3170 
A3180 
A3190 
A3200 
A3210 
A3220 
A3230 
A3240 
A3250 
A3260 
A3270 
A32B0 
A 3290 
A3300 
A 3x3 1 0 
A3320 
A3330 
A3340 
A3350 
A3360 
A3370 
A 3380 
A3390 

IS IN FRROPr NK=rlT4) A3400 

A 34 1 0 
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765 IF(NK~2> 762,766.766 A3420 

766 nCONX=QCON(NK-l>+< ( GcON ( NIK J -OCON ( NK-1 ) ) / ( TI ME ( NK ) -TI MF ( Nk-1 ) ) 1 A3430 

1 * ( T-TTME ( NK-1 ) > A3440 

OC0NX=FC0NV*QC0NX A3450 

(5RADX=QRAn(NK-l > + ( < OR AD ( NK >-OR AD< NK-1 1 ) / ( TIME ( NK )— TIMEINk— 1 ) ) ) A346 0 

1*<T-TTME(NK-1 1 > A3470 

ORADX=FRAn*0RADX A3480 

WFLX=VEL (NK-1 )+( l VEL(NK)-VEL(NK-1 » >/( TIME ( NK >-TImF < NK-1 > ) 1 A3490 

1* ( T — T TMF < NK— 1 ) ) A3500 

60 TO 775 A3510 

770 oCONX=FCONV*OCON(NK) A3520 

ORADX=FRAn*0RAD(NK ) A3530 

VELX=VEL(NK) A3540 

C A3550 

C COMPUTE HFAT RLOCKAGe AT FRONT SURFACE A3560 

775 IF ( I 17—1 ) 77ft, 77ft. 776 A3570 

776 T F ( 1 17— NHP ) 777,777,778 A3580 

777 TF(TX2( INT>— TW( 117) ) 782.788,780 A3590 

778 wF ITE ( 6 , 779 ) TX?(INT) A36nn 

779 FORMAT ( 1H0 , 8i)M THE RANGE OF THE ENTHAL PY-TEMPFR AtURE CURVE FIT WAS A3610 

IeXCFEDED AT A TEMPERATURE 0F,!Eln.4) A36?0 

GO TO 906 A 3630 

780 117=117+1 A3640 

GO TO 77b A3650 

782 IF <TX?< INTI-TW t 117-1 ) ) 7flU,78ft,7ft6 A3660 

784 U7=I17-1 A3670 

GO TO 776 A3680 

786 HW=HX( 117-1 ) + ( ( HX ( 1 1 7 ) — HX (117-1 ) ) / ( TW ( 1 1 7 ) -TW ( 1 17-1 ) ) ) A3690 

1*(TX2(INT )-Tw< 117-1) ) A3700 

GO TO 789 A3710 

788 HW=HX(I17) A3720 

789 hTX=H300+( (VEI X**? ) /60056 • 5 ) A3730 

OPLOCK=(FRLOw*XMOG(lNT)*(HTX-HW) 1/3600. 0 A374P 

C A3750 

C COMPUTE HFAT TN DUE TO SURFACE COMBUSTION A3760 

XMDO=XMDC A3770 

CALL OXIUATIXMDO.QOXTD) A3780 

C A3790 

C COMPUTE Q-HOT WALL A3800 

IF(TDMP.EO.O. 1 GO TO 4001 A3810 

TFIT.GE.TDMP) DMpzl.n A3B20 

4001 7= ( HTV— Hw 1 / ( HTX— H300 1 A3830 

TF(Z-l.O) 790.702,793 A3840 

790 IF ( Z ) 791,791,793 A3850 

791 OHW=0 • 0 A3860 

GO TO 1790 A3870 

792 «HW=OCONX A 3880 

GO TO 1790 A3B90 

793 OHW=Z*QC0MX A3900 

1790 7ZZ= l OHW— OBLOCK ) /QHW A3910 

IF ( ZZ7— 0 . ? ) 1798, 179ft, 1794 A39?0 

1798 oPLOCK=0 ,8*GHW A3930 

C A3940 

C NET HEAT INTO FRONT SURFACE A3950 

1794 IF<IEM(InT)) 795.795,797 A3960 

795 IF(TX8( INTl-TCHAR) 706,796,797 A3970 

796 E m X=EMV A3980 
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GO TO 798 A3990 

797 TFM(lNT)=t A400H 

FMX=EMC A4010 

798 OlN=QRADx+QHw+G0XlD-QRL0CK-(4,8333E-l3)*FMX*FV* ( <TX2( TNT)**4)- A4020 

KTV**U)) A4030 

IF(DMP) 804*804*800 A4040 

800 WRITE ( 6 * 801 ) A4050 

801 FORMAT ( /// ) A406O 

WRITE (6* 802) OCONX , OR ADX , VELX , HT X * HW QRLOCK * QHw * OOX 1 0 * OTN A4070 

802 F0RMAT< 1 X,6HQC0NX=* 1PF1?.5*2X,6HoRAOX=* 1PF12.5* 2X*5HVFLX = , 1PF1 2.5* A4080 

l?y # 4HHTX= t lPFl2.5*2X # 3HHw=*lPF12.5/lX,2H7=* 1PF1?.5*2X*7H0PL0CK=*1P A4090 

2F12.5*2X*UH0HW=* lPEl2.5 f ?X*6H00XiD=*lPEl?.5*2Xr4H0IN= f 1PF12.5/) A4100 

804 OTN=Q!N*3600,0 A4110 

A4120 

CHECK FOR PROMT SURFACE RFCESSION (CHAR LAYER RFmOVAL) A4130 

CALL RECESS ( XWDC *XLOST*TRFC*DT*PHOC*TS* SR *TX2( 1 ) , NRFC * NRS * ERRS * Sy 0 A4140 

1 « SDOT * DMP ) A 4 150 

TF(ERRb) 8050*8080*9(18 A4160 

8050 VLV=VLV~XI OST A4170 

XI STV=XLSTV+Xl OST A4180 

XI STI=XLSTV*12.0 A4190 

nxyrVLV/ ( XNP-t ,0) A4200 

XV< 1 ) =0 * 0 A4210 

no 1780 I=2*NP A4220 

XV < D=XV( T-l )+OXV A 4230 

1780 C.ONTIMUF A4240 

PX=DXV A4250 

T F ( FRR4 ) 806,806.805 A4260 

80b G 0 TO 90b A4270 

806 CALL COFFF ( NpFT * SHOT ) A4280 

TF(DMP) 8069,8069*8081 A4290 

8061 wP I TE ( 6 * 8062 ) A4300 

8062 FORMAT ( / 1 X * 23H COFFFICIFNTS FOR SWUFT/) A4310 

nO 8066 X=1*NPFT A4320 

WRITE (6* 8064) A(I)*R(I),r(I)*D(I)*I A4330 

8064 FORMAT U HO *5HA ( I ) = * 1PF12 . 5 * 2X *5HR ( I ) = , 1PF1 2 . 5 * 2X , 5HC ( T ) = • 1 PF1 2 • 5 . 2 A434 0 

1X*5HD( I )=*1PF1 2.5*2X, ?HI=*I3) A4350 

8066 CONTINUE A4360 

8069 T F ( FRR2 ) 807*807,805 A4370 

807 TFIFRR3) 810,810*808 A4380 

808 WRITE ( 6 * 809 ) IKK A4390 

809 FORMATMHO*l8M THF VaLUF OF IKK = ,1I4) A4400 

GO TO 90b A4410 

810 CALL SWUFT (A,R*C*D*Ty* NPFT * DMP ) A4420 

827 no 828 1=1 *NP A4430 

TX 1 ( I ) =TX2 ( I ) A4440 

TX2(I)=TY(I) A4450 

828 CONTINUE A4460 

CALL nON2 ( XLOST *XV*TX2* NP * XC , TX2C * XDV , KK \/ * XLST V • HXX ) A4470 

830 CALL ABLATE A4480 

XMDT=XM0G( INT)+XMnC A4490 

LT=NP+1 A4500 

no 1815 1=1*NMB A4S10 

LLT=NPM(1) A4520 

TF(I.FQ.l) GO TO 1812 A4530 

XF(GAPX( X-l) .FQ.O* ) GO TO 1812 A4540 

KKT = 1 A 4550 
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GO TO 1813 A4560 

1812 KKT=2 A4570 

1813 no 1815 J=KKT,LLT A4580 

TX2T( Jr I )=TY(LT) A4590 

l T=LT+1 A4600 

1815 CONTINUE A4610 

no 1819 irlrNMB A4620 

TF(I.FQ.l) GO TO 1 81 8 A4630 

TF ( GAPX { 1-1) . FQ • 0 • ) GO TO 1817 A4640 

GO TO 1810 A465P 

1816 TX2TM rl )=TY(WP) A4660 

GO TO 1819 A4670 

1817 IX=NPM(T-1) A468D 

TX2T l 1 r I )=TX2T (LXr 1-1 ) A4690 

1819 CONTINUE A4700 

I M=N P + 1 A4710 

no 833 T=1 A4720 
L7=NPM(I) A473H 

nO 833 J=1 rL? A4740 

TX2(LM)=TX2T( Jr I ) A4758 

l M=LM+1 A4760 

833 CONTINUE A4770 

00 5834 I -2 r NpTT A4780 

TP ( T -TT ABLE ( I ) ) 5835,B835r5834 A479D 

5835 nTS=DFUTT ( 1-11 A4800 

TPRCT=IPKC (1-1 ) A4810 

OT=OELTT ( T-l ) / 360 0.0 A4820 

GO TO 5836 A4830 

5834 CONTINUE A4840 

nTS=UELTl (NPTT) A4850 

TPRCT=IPKC (NPTT) A4860 

OT = DEL.TT (NPTT) /360 0, o A4870 

5836 TCT=ICT+1 A4880 

5838 Vl.TEM=SAVY3 A4890 

CALL TS0THM(XVrTX2r 1060. r NP r SAVE I T ) A4900 

FAvEITrSAVElT+XLSTV A4910 

TF (SAVY3.LT .SAVFIT) F a VY3=FAVF TT A4920 

TF(VLTEM # FQ.SAVY3)G0 TO 839 A4«30 

FAVX=T A4940 

FA VY1-XLST I A4950 

FAVY2=TV2 (NP) A496D 

CALL TSOI HM(XVr TX?r 1460. .MPrSAVYu) A4970 

839 R| TFM=SAVY4X A498D 

CALL TSOI HM( XVr TX2r 1 46 0. r NP r WFKFpP ) A4990 

wFKFEP=WEKEFP+XLSTV A50n0 

TF( SAVY4X.LT. WEKEFP)c;AVYuX=WFKEFP A50) n 

TF(BLTEM.F0.SAVY4X)G0 TO 838 A50?n 

fAVEXV^T A 5 0 3 0 

FA VY 1X-XLST I A504D 

FAVY2X=TX? (NP) A5050 

CALL TSOTHM (XVrTX2tlfl60. * NP r F A VY3X ) A 5060 

838 CONTINUE A5070 

TF(IPRCT-TCT) 835 » 835 r 84n A508D 

835 WRITE (6r 837) T r OCONX , OR A OY , VFt X r XMDG ( TNT) r XMDO r X mHT r XL ST I r QHW A5090 

837 FORMAT ( 1H0 r 5HT I WE- r A5lnn 

1 1PF1 2.5r2X, 12HQC0NVFCTTVE=r 1PF1 2.5. ?Xr 1 1H0RADTAT A51 10 

lTVE=r lPF12.5r5>Xr9HVF L 0CTTY=rlPE12.5/lXr 18HGAS ARlATION RATEr.lPFi? A5120 
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2 . 5 * ?X , 19HCHAR ABLATION R ATE= , 1 PFl 2 . 5 . ?X , POHTOT AL ABLATION RATF=,lP A5130 
3F!2.5/1X,16 HrFCFSSI0n DFPTH=,1PFl?.5,?X,10HOHOT wALL=,1PE12.5) A5140 

840 T-T+DTS A5150 

841 IFCNPLOT.NE. 1 ) GO TO 842 A5160 

CALL SAVE( ASAVE1, ASAv/F2.» AS AVF3 * Uc;F A , NX A , XLST I , DT5 , TL I M . T , V ALUFA J A5170 

CALL SAVE<BSAVEl,RSAvF2,BSAVF3.U5FB»NXB,TX2(NP) ,DTS»Tl IM,T,VALUFR) A5180 
CALL TSOTHM(XV*TX?» 1060, ,NP,Y3) A5190 

CALL RAVE ( C5 AVE1 , CS AvF2 * CS AVF3 , USF C,NXC,Y3* DTS » TL IM > T * V ALUEC ) A5200 

CALL IS0THM(XV.TX2* 1460. ,NP.Y4) A5210 

CALL SAVE(DSAVE1 ,DSAvF2»DSAVF3, USFD,NXD,Y4 , DTS. TlIM,T, VALUED) A5220 

T F ( USF A . NF .0.0) GO TO 9842 A5230 

I F ( USFB . NF . 0 . 0 ) GO TO 984 2 A5240 

T F ( USFC . NF .0.0) 60 TO 984? A5250 

T F ( USFO • NF . 0 . D ) GO TO 9842 A5260 

GO TO 9843 A5?70 

9842 xPLOT=T-DTS A5280 

Y PLOT 1 -VALUE A A5290 

iF(usfa,nf.o.o)Yploti=u9fa assoo 

YPL0T?=VALUEb A5310 

TF(USFB.NF.O.O) YPLOT?=USFP A53?0 

YPL0T3=VALUFC A 5 3 3 n 

T F ( USFC . NF . 0 . 0 ) YPLOT 3=USEC A5340 

YPLOT AzVALUFn A5350 

TF(USFD.NF.O.O) YPL0T4=USFH A5360 

WRITE (11)XPL0T,YPL0T1 ,YPI 0T2,YP| 0T3,YPLOT4 A5370 

9843 TF(ICTP.NF.O) GO TO 842 A5380 

T C TP— 1 A5390 

XPLOT=T A5400 

YPL0T1=XLSTI A5410 

YPLOT?=TX?(NP) A5420 

CALL TSOTHM (XV,TX?»lO60, , NP , YPL0T3 ) A5430 

CALL TSOTNMI XV,TX?» 1460, ,NP, YPL0T4) A5440 

WRITE (11)XPloT»YPLOt1 * YP1.0T2 , YPL0T3, YPL0T4 A5450 

842 rF(IPPCT-TCT) 845,845,900 A5460 

845 WRITE ( 6 * 850 ) T A5470 

IPCT=TPCT+1 A5480 

TFIIPCT.E0.2) TPCT=0 A5490 

TF(IPCT.EO.O) TCTPrO A5500 

850 FORMAT ( 1 HO » 7PHTFMPER ATURF DISTRIBUTION IN HFAT SHIELD AT THF FND 0 A5510 

IF THE TIMF STFP, T= ,1 PE)?. 5, IX, 7HSFC0NDS// ) A5520 

wP ITE ( 6 , 860 ) A5530 

860 F0RMAT(4X,49HTEMPFRATURF DISTRIBUTION IN THF ABLATING MATFRI Al // t A5540 

kKV=KKV+1 A5550 

wP I TE ( 6 * 862 ) (TX?C ( I ) , 1 = 1 ,KKV) A5560 

862 FORMAT ( 6X » 1PF 12.5, 1P5F 16 . 5 ) A5570 

TJ=NP+1 A5580 

wR I TE ( 6 , 864 ) A5590 

864 F0RMAT(//4X, aoHTRMPFr ATURF DISTRIBUTION IN THF BACK-UP STRUCTURF// A5600 

1) A 56 10 

wR ITE ( 6 * 862 ) (TX2 t I ) , T = IJ,NPF) A5620 

WRITE 1 6 r 865 ) A5630 

865 FORMAT ( // ) A5640 

TCT=0 A5650 

900 CONTINUE A5660 

TF(T-TLIM) 750*750,905 A5670 

905 IF(NPLOT.NE.l) GO TO 909 A5680 

X AVY3=SAVY3-SA VYl/12 , A5690 
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XAVY4X=SAVY4X-SAVY1X/12. A5700 

TFISAVX.EQ.XPLOTIGO TO 9005 A5710 

wPITE ( 1 1 ) SAVX , SAVY1 . SA VY2 . X AVY3 . S AVY4 A5720 

9005 TFlSAVEXX.EQ.XPLOTIGO TO 9006 A5730 

SAV4I=SAVY4X*l2. A5760 

9006 SAV3I=SAVY3*15>. A5750 

WRITE ( il > SAVEXX*SAVYlX* SAVY2X .SAVY3X* XAVY4X A5740 

wPITE(6.929>SAV3I.SAv4l A5770 

9?9 FORMAT 1 1H0 * 23HMAXIMUM 1060 ISOTHFRM =E16.8*?X23HmAXIMuM 1460 ISOTH A5780 
1FPM SF16.R) A5790 

WRITE ( 1 1 ) STOP * STOP * STOP . STOP * ST OP A5800 

909 TP (LPLOT .NE.NCASEIGO TO 911 A5810 

DATA FND/6H FND / A5820 

wP I TE ( 1 1 > FND , FND . FND . FNn . FND . FND, FND* FND . FND . FND. FNO • FND A5830 

OllIT=888tt. A5840 

WRITE ( 1 1 )QUIT* QUIT* QUIT* QUIT .QUIT A5B50 

FND FILE 11 A5860 

RFWIND 11 A5R70 

911 IF ( TEST2 ) 910.930.930 A5880 

910 DO 920 JUK=1,NPF A5890 

TX1 ( JJK)=TEMDT ( JJK ) A5900 

TX2<vJJK)=TX1(JJK) A 59.10 

TUL1(K)=TX1(K) A5920 

Tl)L2(K)=TXl(K) A5930 

9?0 CONTINUE A5940 

IL=NP+1 A5950 

no 926 1=1* NMR A5960 

Tl N=NPM< I ) A5970 

00 924 J=l*Il_N A5980 

TX 2 T(J.n=TFMDl(IL) A5990 

Tl=IL+l A6000 

9?4 CONTINUE A6010 

9?6 CONTINUE A6020 

BO TO 940 A6030 

930 CALL TEMPO A6040 

940 T=TINT A6050 

nx=0xx A6060 

nTS=DFLTT(l) A6070 

nTsDELTT ( 1 > /.3600 . 0 A6080 

VI. V=VL A6090 

GO TO 50 A6100 

END A61 tO 
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SIBFTC COEF 

C THIS SUBROUTINE DETERMINES THE COEFFICIENTS OF THE MATRIX 

SUBROUTINE COEFF (NPFT , SDOT > 

C 

DIMENSION TITLE (12) *HEADNG (12) r XIDNT (12*12) * TKC ( 20 ) * XKC ( 20 ) » 
1CPC(20) * TKV (20 ) *XKV(20) * TCPV (20 ) * CPV ( 20 ) .TIME ( 300 ) * QCON ( 300 ) , 

2QRAD ( 300 ) *VEL(300) * XNPM ( 12 ) * NKPB ( 12 ) *NCPB(12) *TXK(20,12) *XK(20,12) 
3 * TCP (20 * 12) * CPX ( 20 * 12 ) ,RHOBX(12) *XBM(12) »EMFB(12) *EMBB(12) *HXX(12) 
4 * GAPX ( 12 ) *FTEST ( 12 ) *BTEST(12) * TEMDI ( 200 >* TX1 (200 ) *TX2(200) * 

5TX2T (10*12) * TUL1 (200 ) * TUL2 (200 ) *HX ( 50 ) * Tw ( 50 ) . IR (50 ) , IR1 ( 50 ) , 

6IR2 ( 50 ) * TUL ( 50 ) * I EM ( 50 ) * TY ( 200 ) * A (200 > * B (200 ) * C (200 ) , D < 200 ) * 

7R ( 50 ) * RHO ( 50 ) *CP(50) ,DXB<12) *XKB(10*12) *CPB ( 10 * 12) * XMDG ( 50 ) * 
8YK(50) * AB (10*12) *BB ( 10 * 12 ) *CB ( 10 . 12) *DB(10*12) *SB(10*12) * 

9RB1 (10*12) *RB2(10*12) *H ( 12 ) * S ( 50 ) * NPM ( 12) 

DIMENSION TTUL ( 50 ) * RHOY1 ( 50 ) * RH0Y2 ( 50 > , DRHO ( 50 ) * TCPC ( 20 > 

C 

COMMON TKC * XKC * TCPC * CPC * TKV * XKV * TCPV * CPV . XNPM * RHOBX r XBM * EMBB * 

1EMFB * NKPB* NCPB * TXK *XK* TCP * CPX * NPM * GAPX * FTEST *BTEST * TEMDI *TX1* 
2TX2*TX2T * TUL * TUL1 * TUL2 *IR*IR1*IR2*A*B*C*D*S*R*AB*BB*CB*DB*SB* 

3RB1 *RB2*TY *RH0Y1 * RH0Y2 * XMDG *RHO *CP * YK * XKB*CPB *DXB * DT * XLOST * 

4TABL * TCHAR * TREC * RHOV * RHOC * FBLOW * EMV * EMC r H300 * NKC * NCPC * NK V * NCPV * 

5NP * NMB * NPBS » NPF * TEST2 * TEMPI *TX0 * TENV * HENV * FENV * QLOSS* TLIM * TINT 
COMMON Ilf 12*13 * 14 * 15* 16 * QIN* INT *DX * XMT * TL* VL *BL *DMP* ERR1 * ERR2 * 
1ERR3*ERR4* HV* VPT* CHARKrCHARC* ABLK* ABLC* XMDCrH 
C 

CALL PROP 
YNP=NP 

S(INT)=(RHO( INT)*DX*CP(INT> )/(2.0*DT> 

R(INT)r(i.0)/( ( DX/2 . 0 ) * ( ( 1 . 0/YK ( INT ))+( 1 . 0/YK ( INT+1 > ) ) > 

A( INT)=0,0 

B ( INT ) = ( - ( ( XMDG ( INT ) +XMDC ) *CP ( INT ) +S ( INT ) +R ( INT ) -RHO ( INT ) *CP ( INT ) 
1* (SDOT/ ( 2. 0*( YNP-1. 0) ) ) ) ) 

C ( INT ) =XMDG ( I NT+1 )*CP( INT+1 ) +R ( INT) +RHO( INT+1 )*CP( INT+1 )*SDOT 
1*( ( YNP-1 , 5 ) / ( YNP-1 • 0 ) ) 

D(INT)=(-(QIN+S(INT)*TX2(INT) ) ) +( XMDG ( INT) -XMDG ( INT+1 )) *HV 

NPP=NP-1 

JNT=INT+1 

DO 10 I=JNT»NPP 


XI = I 
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S ( I ) = ( RHO ( I ) *DX*CP ( I ) ) /DT 

R(I)=(1.0)/( (DX/(2.0*YK(I> ) )+(DX/(2.0*YK(I+l) ) ) ) 
A(I)=R(I-1) 

B ( I ) = ( - ( XMDG ( I ) *CP ( I ) +R ( 1-1 ) +R ( I ) +S( I ) + RHO ( I ) *CP ( I 
1-0 . 5) / ( YNP-1 • 0 ) ) ) ) 

C ( I ) =XMDGjl+l-)-*CP_a + l ) +R ( I ) +RHO ( 1 + 1 ) *CP ( 1 + 1 ) +SDOT* ( 

1/ ( YNP-lriOJ-J — - ~ ^ ' — ^ 

D ( I ) = (-<S CD *TX2<T) ) ) +f XMDG ( I) -XMDG ( 1 + 7) ) *HV 
CONTINUE 

R(NP)=(1.0)/( (DXB(l) / ( 2 • 0*XKB (1*1) ) ) + (DXB ( 1 ) /( 

S ( NP ) = ( RHO ( NP ) *DX*CP ( NP ) +RHOBX ( 1 ) * CPB (1*1) *DXB ( 1 ) ) 
A(NP)=R(NP“1) 

B(NP)=(-(XMDG(NP)*CP(NP)+R(NP-1)+R(NP)+S(NP) ) ) 
C(NP)=R(NP) 

D(NP)=(-S(NP)*TX2(NP) ) +XMDG(NP)*HV 
DO 200 I=1*NMB 
IF(I-l) 20*20*30 
AB(1*I)=A(NP) 

BB ( 1 * I ) =B ( NP) 

CB(1*I)=C(NP) 

DB(1*I)=D(NP) 


) *SD0T* ( ( YNP-XI 
(YNP-XI-0.5) 

2 • 0*XKB (2 * 1 ) ) ) ) 
/ ( 2. 0*DT > 


B0000 

B0010 

B0020 

B0030 

B0040 

B0050 

B0060 

B0070 

B0080 

B0090 

B0100 

B0110 

B0120 

B0130 

B0140 

B0150 

B0160 

B0170 

B0180 

B0190 

B0200 

B0210 

B0220 

B0-230 

B0240 

B0250 

B0260 

B0270 

B0280 

B0290 

B0300 

B0310 

B0320 

B0330 

B0340 

B0350 

B0360 

B0370 

B0380 

B0390 

B0400 

B0410 

B0420 

B0430 

B0440 

B0450 

B0460 

B0470 

B0480 

B0490 

B0500 

B0510 

B0520 

B0530 

B0540 

B0550 

B0560 

B0570 

B0580 

B0590 
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GO TO 65 B0600 

30 L=NPM(I-1) B0610 

IF (FTEST ( I ) ) 45»40»45 B0620 

40 SB(l»I)=(RHOBX(I>*CPB(l»I>*DXB(I>+RHOBX(I“l)*CPB(L»I-l)*DXB(I-l> )/ B0630 

1(2.0*DT) B0640 

RB1(1»I)=(1.0)/((DXB(I-1)/(2.0*XKB(L#I-1)))+(DXB(I-1)/(2.0*XKB(L-1 B0650 

1»I“1)))) B0660 

RB2(l»I)=(1.0)/( (DXB(I)/(2.0*XKB(1>I) ) ) + (DXB ( I ) / (2 . 0*XKB (2 » I ) ) ) ) B0670 

AB ( 1 r I ) =RB1 ( 1 » I ) B0680 

BB(1»I)=(“( RBI ( 1 rl ) +RB2 ( 1 » I ) +SB (1*1)) ) B0690 

CB(1*I)=RB2(1»I) B0700 

DB(1»I)=(-(SB(1>I)*TX2T(1»I) ) ) B0710 

60 TO 65 B0720 

45 IF (FTEST ( I ) ) 50»40*55 B0730 

50 G=(1.73E-09)/(1.0/EMBB(I-1>+1.0/EMFB(I)“1.0) B0740 

60 TO 60 B0750 

55 6=0 • 0 B0760 

60 SB ( 1 r I ) = (RHOBX ( I ) *CPB ( 1 » I ) *DXB ( I ) )/(2.0*DT) B0770 

RB2(l,I)=(1.0)/( (DXB(I)/(2.0*XKB(1»I> ) ) + (DXB ( I ) / (2. 0*XKB (2 » I > ) ) ) B0780 

AB(lrI)=H(I-l)+4.0*G*(TX2T(LrI-l)**3) B0790 

BB ( 1 ► I ) = ( - ( H ( 1-1 ) +4 . 0*G* ( TX2T ( 1 » I ) **3 > +RB2 ( 1 1 1 ) +SB ( 1 > I ) ) ) B0800 

CB ( 1 r I ) =RB2 ( 1 r I > B08X0 

DB(lrI)=3.0*G*( (TX2T(L» 1-1 > **4 ) - (TX2T ( 1 » I > **4 > ) -SB ( 1 » I > *TX2T ( 1 , I ) B0820 

65 LF=NPM(I)-1 B0830 

DO 100 J=2»LF B0840 

SB (J r I ) = (RHOBX (I)*CPB(J»I) *DXB ( I ) ) /DT B0850 

RBl(J»I) = (1.0)/( (DXB(I)/(2 • 0*XKB ( J“1 r I ) ) ) + (DXB ( I ) / (2. 0*XKB ( Jr I ) ) ) ) B0860 

RB2(J,I)=(1.0>/( (DXB(I)/(2.0*XKB(J+1»I) ) > + (DXB ( I > / (2. 0*XKB ( J» I ) ) ) ) B0870 

AB(JrI)=RBl(J'I> B0880 

B8(J»I)=(-(RB1(J>I>+RB2(J»I)+SB(J>I) ) ) B0890 

CB < J » I ) =RB2 ( J r I ) B0900 

DB( J»I ) = (-(SB( I> *TX2T( Jr I) ) ) B0910 

100 CONTINUE §0920 

IF(I-NMB) 110 » 250 > 250 B0930 

110 LNF=NPM(I) B0940 

IF (BTEST ( I ) ) 120 r 115* 120 B0950 

115 SB ( LNF • I ) = ( RHOBX ( I ) *CPB ( LNF 1 1 ) *DXB ( I ) +RHOBX ( 1 + 1 ) *CPB ( 1 1 1 + 1 ) *DXB ( 1+ B0960 

11) )/(2.0*DT) B0970 

RBI (LNF *I)=(1.0)/( (DXB(I)/ ( 2. 0*XKB (LNF-1 * I ) ) ) + (DXB ( I ) / (2. 0*XKB (LNF B0980 

X r I > > ) > B0990 

RB2(LNF*I)=(1.0)/( (DXB ( 1 + 1 ) / (2. 0*XKB ( 1 * 1 + 1 ) ) >+(DXB(I+l)/ B1000 

1(2.0*XKB(2*I+1) ) ) ) B1010 

AB (LNF» I ) =RB1 (LNF • I ) B1020 

BB (LNF* I ) = (- (RBI (LNF * 1 ) +RB2 (LNF * 1 ) +SB (LNF * I ) ) ) B1030 

CB ( LNF * 1 ) =RB2 ( LNF * 1 ) B1040 

DB(LNF*I)=(-(SB(LNF*I)*TX2T(LNF*I) ) ) B1050 

60 TO 200 B1060 

120 IF (BTEST ( I ) ) 125* 115* 127 B1070 

125 G=(1.73E-09)/(1.0/EMBB(I>+1.0/EMFB(I+l>-1.0) B1080 

60 TO 130 B1090 

127 6=0.0 B1100 

130 SB (LNF * I ) = (RHOBX ( I ) *CPB (LNF * I ) *DXB ( I ) ) / (2* 0*DT ) B1110 

RBI (LNF *I)=(l«0)/( (DXB (I) / (2. 0*XKB (I .NF-1 * I ) ) ) + (DXB ( I ) / ( 2. 0*XKB (LNF B1120 
lrl)))) B1130 

AB(LNF * I ) =RB1 (LNF * 1 ) B1140 

BB (LNF ► I ) = (- (RBI (LNF * I ) +H ( I ) +SB (LNF * I ) +4 . 0*G* ( TX2T (LNF ► I ) **3) ) ) B1150 

CB(LNF r I )=H( I ) +4. 0*G* (TX2T (1*1+1) **3) B1160 

DB(LNF*I)=3.0*6*( < TX2T ( 1 * 1 + 1 ) **4) - ( TX2T (LNF * I ) **4) ) “SB (LNF ► 1 ) *TX2T B1170 

1 (LNF * I ) B1180 

200 CONTINUE B1190 
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250 MN=NPM(NMB> B1200 

IF (QLOSS) 270 » 260 , 270 B1210 

260 SB t MN , NMB ) = ( RHOBX (NMB ) *CPB (MN ,NMB) *DXE ( MMB ) )/(2.0*DT) B1220 

KBl(MN,NMB)=(1.0)/( ( DXB (NMB ) / (2 . *XKB (MN-1 , NMB ) ) ) + ( DXB ( NMB ) / ( 2. 0*XK B1230 
IB ( MN , NMB ) ) ) ) B1240 

AB(MN,NMB)=RB1 (MN»NMB) B1250 

BB (MN » NMB) = ( - (RBI ( MN , NMB) +SB (MN , NMB) ) ) B1260 

CB(MN,NM6)=0.0 B1270 

DB (MN , NMB ) = ( - ( SB ( MN ,NMB ) *TX2T (MN, NMB) ) ) B1280 

GO TO 280 B1290 

270 SB ( MN » NMB) = ( RHOBX ( NMB ) *CPB ( MN , NMB ) *DXB ( NMB > > / ( 2. 0*DT ) B1300 

RBI ( MN , NMB > = ( 1 • 0 ) / ( ( DXB < NMB ) / ( 2 . 0*XKB ( MN-1 , NMB ) ) ) + ( DXB ( NMB) / ( 2 . 0*X B1310 

1KB (MN» NMB) ) ) ) B1320 

AB(MN#NMd)=RBl(MN»NMB) B1330 

BB ( MN , NMB ) = ( - ( RBI ( MN > NMB ) +HENV+ < 1 . 73E-09 ) *FENV*4 . 0* ( TX2T ( MN » NMB ) ** B1340 

13 ) +SB ( MN » NMB ) ) ) B1350 

CB ( MN t NMB ) = 0 • 0 B1360 

DB(MI\HNMb)=<-(HENV*TENV+FENV*<1.73E-09)*< (TENV**4)+3.0*(TX2T(MN»NM B1370 

IB ) **4 ) )+SB(MN»NMB)*TX2T(MN»NMB) ) ) B1380 

260 L=NP+1 B1390 

DO 300 1=1 » NMB B1400 

K=NPM ( I ) B1410 

IFCI.EQ.l) GO TO 282 B1420 

XF(GAPXd-l) .EQ.O.) GO TO 282 B1430 

KT=1 B1440 

GO TO 28S B1450 

262 KT=2 B1460 

265 DO 290 J=KT,K B1470 

A(L)=AB(o*I) B1480 

B(L) =BB ( J f I ) B1490 

C(L)=CR(Jd) B1500 

D ( L ) =DB ( J f I ) B1510 

IF(DMP) 289,289.286 B1520 

28b WR ITE ( 6 ,287 ) AB ( J, I ) , BB ( J , I ) , CB ( J > I ) , DB ( J, I ) , J, I , A (L ) , B ( L ) , C (L) , D ( B1530 

1L),L B1540 

287 FORMAT (1H0,8HAB( J, I > = , 1PE12. 5 » 2X » 8HBB < J, I ) = , 1PE12 . 5 , 2X , 8HCB ( J . I)=. B1550 

1 1PE12 .5 » 2X» 8HDB ( J , I ) = , 1PE12 • 5 , 2X , 2HJ= , 13, 2X , 2HI= , 1 3/1X , 5HA ( L ) = , 1PE B1560 

212.5,2X,5HB(L) = * 1PE12 . 5 * 2X * 5HC (L ) = , 1PE12 . 5, 2X , 5HD (L ) = , 1PE12 . 5 . 2X . 2 B1570 

3HL= ,13) B1580 

289 L=L+1 B1590 

290 CONTINUE B1600 

300 CONTINUE B1610 

NPFT=L— 1 B1620 

RETURN B1630 

END B1640 
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IbFTC PRP 

THIS SUBROUTINE DFTFrminfs THF PHYSICAL PROPERTIfS OF THf 
HEAT SHIELD STRUCTURE 
SUBROUTINE PROP 

DIMENSION TITI E(l?) r HFAPNO < 1? ) # X T DNT ( 1 2 , 1 ? ) #TKC(?0) #XKC(20) , 
lCPC(20) #TKV(20) rXKV(pn) ,TCPV(20 ) , CPV ( 20 ) , TI ME ( 30f) ) #GCPN<300) , 

20PAD ( 300 ) tVEL (300 ) # XnPM < 1 ? > » NXPB ( 1 2 ) ,NCPR( 12) #TXk <20, 12) ,XK<20# 1?) 
3, TCP<?0#12) #CPX(?0# t?) # RHOBX (12) , VBM (12) ,FMFR( 12) #EMrr< 1?) ,HXX< 1?) 
4 , GAPX ( 12 ) ,FTFST(1?) # flTEST (12 ) »TFmPI ( 20 0) , TX 1 (200) #TX?(20n) # 

5TX2T( 1 Of 1?) #TUL1 (200) ,TUL?(200) # HX < 50 ) • Tw < 50 ) • IR ( 50 > , TR1 < 50 ) , 

6 T P2 ( 50 ) # ruL(SO) # iEM(sn) ,TY(200) f A ( 200 ) #R(?00> #C<200) ,0(200) , 

7R ( SO ) #RHO (50 ) ,CP<50> ,OXP< 12) #XKB<10#12) #CPB(10#l2) # XMOG(SO) , 

BYK (50) # Ah (10. 12) ,RB< 1 0, 12) ,CR(10, 12) # OB < 10» 1?) # SR (10, 12) . 

9rP1 <10#12) #RR?< 10# 1?) #H( 1?) #S(50) #NPM(1?> 

DIMENSION TTuL ( 50 ) #RH0Y1 (SO ) .RHOY?(50 > #DrHO(50) #TCPC(?0) 

common TKC# XKC # TCPC#CPC,TKV#XKV#TCPV,CPV, YNPM,RHoPX,XPM#FMBB, 
lFMFRfNKPBf NCPR#TXK #XK ,TCP,CPX,NPM#GAPX,FTEST,PTEsTf TFMDI,TXI , 
2TX?,TX2T,TUL,TUL1 #TUL?# IR#IR1 # IR? # A # B , C # D # S # R , AP , RB # CP # DH # SP , 

3RP1 #RR2# TY,RhOY 1 ,RHOy?# XMPG# RHO, CP# YK,XK r,CPB,DXB#DT,yLOST, 

4tabl,tchar#trfc#rhov,rhoc#fplow#fmv#emc#h30o#nkc,ncpc,nkv#ncpv# 

5NP • NMR # NPRS # NPF , TEST 2 # TEMP I # TXO # TENV # HENv , FFNV , Ol.OSS , TLI M • T INT 
COMMON II #12. 13# I4#I5. IE,OlN,TNT,nX#XMT#TL#VL#BL,PMP,FRPt #ERR?, 
1 fRR3#FRR4#HV, VPT.CHARK.CHARC# ABLk# ABLC#XMPC#H 

KTNT=TNT 

no 170 i=kint#np 
)0 I F ( IR ( I ) ) 12 , 1 2 # 1 00 
12 TUL ( I )=AMAX1 (TXK I ) #TX2( T ) ) 

TF(TUL(I) .LE.TABL) go TO 20 
T P ( I ) -1 
GO TO 10o 

20 T F ( 1 1—1 ?5 # ? #2i 

21 TF(Il-NKV) ?2,a?,P 5 

22 TF(TX?(I )-TKV< ID ) 3S#55,30 

25 WRITE ( 6 # 26 ) TV2(I) 

26 FORMAT ( 1 HO # 87H THE RaNGF OF ONE OF THF ABLATION PROPERTY CURVF FTT 
IS WAS EXCFEDEP AT A TFMPFRATUPE OF #1PE1?*5) 

FRR2=1.0 
GO TO 355 
30 Tl=Il+l 
GO TO 21 

35 TF(TX?( I )-TKV< 1 1 — 1 ) ) 40 # 55 # 50 
40 T 1 =11-1 
GO TO 20 

50 YK (1)=XKV(I1-T ) + < (XKV( I) >-XKV(II-l ) >/(TKV(Il )-TKV(Il-1 ) ) ) 
l*(TX2(I>-TKV(Tl-l) ) 

GO TO 60 
55 YK ( I )=XKV (ID 

60 TF(I2-1> 25 # ?5 # 6 1 

61 I F ( 1 2— NCPV ) 62 # 62 #25 

62 TF(TX2(I )-TCP\M 12) ) 70#P5#65 
65 T2=I2+1 

GO TO 61 

70 TF(TX?<n-TCPV<l2-l) ) 75.85#80 
75 T 2=12-1 


cooon 

room 

C 0020 

C0030 

C0040 

C0050 

CO 060 

CO 070 

COOflO 

C0090 

coion 

c o i in 

CO 120 

CO 1 30 

CO 140 

CO 1 50 

CO 160 

CO 170 

CO 1 80 

CO 190 

C0200 

C021 0 

C 0220 

C0230 

C0240 

C0250 

C0260 

C0270 

C0280 

C0290 

C0300 

C0310 

C0320 

C0330 

C034O 

C0350 

C0360 

C0370 

C0380 

C0390 

C0400 

C0410 

C 0420 

C 0 4 3 0 

C0440 

C045O 

C0460 

C0470 

C0480 

C0490 

C0500 

C0510 

C0520 

C 0530 

C0540 

C055O 

C0560 
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RO TO 60 00570 

AO CP<I>=CPV(I2-1 > + < (CPV(I2)-CPV(I2-1> >/(TCPV<I 2 )-TcPV(I 2 -l) ) ) C05B0 

1*(TX2(I)-TCPV(I2-1> > 00590 

RO TO 90 00600 

A5 CP(I)=CPV(I2> C061O 

90 RHO(I)=RHOV C0620 

RO TO 170 00630 

100 TIJL( I )=AMAX1 <TUL( I > ,TX2( I > > 00640 

TF<TUL(I)-TChAR) 110,110.115 C0650 

110 RHO< I )=RHOV+( RHOV-RHoC >*< <TUL(I>-TARL)/(TABL-TCHAR> ) 00660 

VK ( I ) =CHARK+ ( ABLK— CHARK) * ( ( RHO ( I ) -RHOC ) / ( RHOV-RHOC > ) 00670 

CP ( I ) =CHARC+ ( ABLC-CHARC > * < < RHO < I ) -RHOC ) / ( RHOV-RHOC ) ) 006AO 

RO TO 170 C069O 

115 IF(VPT) 116,116,117 C0700 

116 TTUL ( I > = Tl)L ( I ) C0710 

RO TO 120 00720 

117 TTuL(T)=TX2<I) 00730 

120 TP ( 13-1 ) 25,25*121 C07UO 

121 1F(I3-NKC) 122,122,25 00750 

122 IF < TTMH I )— TkO ( 13) ) 124,135,123 00760 

123 13=13+1 0.0770 

RO TO 121 00760 

124 TFITTULU)— TKC(I3-1> > 125,135,130 00790 

125 T 3=13-1 00600 

RO TO 12 1) C0610 

130 YK ( I)=XKC( 13-1 )+< (XKC( I3)-XKC(I3_1 ) ) / ( TKC ( 13) -TKc ( 13-1 ) ) > 00620 

1 * ( TTUL < I ) — TKC ( 1 3—1 ) ) 00630 

RO TO 140 00640 

135 YK(I)=XKC(I3> 00650 

140 IF ( I 4—1 ) 25,25,141 00660 

141 T F ( 1 4— NCPC ) 142,142,25 00670 

142 TF(TTIIL(I)-TCPC(I4) ) 150,165.145 00660 

145 14=14+1 00690 

RO TO 141 C090O 

160 TF(TTUL(I)-TCPC(I4-1) ) 155,165,160 00910 

155 14=14-1 00920 

RO TO 140 00930 

160 cP ( I ) =CPC (14—1 )+( (CPC( 14) -CPC (14-1) ) / ( TCPC ( 14 ) -TcPC t 14-1 ) ) ) 00940 

1* ( TTUL ( I ) -TCPC ( 14-1 ) ) 00950 

RO TO 166- 00960 

165 CP< I )=CPC ( 14) 00970 

166 RHO(I)=RHOC 00960 

170 CONTINUE C0990 

01000 

DETERMINATION OF PROPFR BACK-UP 5HIFLD MATERIAL PROPERTY 01010 

01020 

DO 300 1=1, NMR 01030 

DXB ( I ) =XBM ( I ) /( ( XNPM ( T 1 — 1 .0)612,0) 01040 

LKP=NKPB(I) 01050 

LCP=NCPB(I) C1060 

NN=NPM(I) 01070 

DO 280 J=1»Nn 01060 

200 IF ( 15-1 ) 2 03,203.201 01090 

201 T F ( 1 5— LKP ) 202,202,203 01100 

202 IF(TX2T( J, I)-TXK( 15,1) ) 206,220,205 OHIO 

203 WRITE (6,204 ) T,TX2T(J,I) 01120 

204 FORMAT* 1 HO, 32H THE RANGE OF ONE OF THE Nl iMBFR ,I2«71H BACKUP FTPllC C1130 
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1TURF PROPFRTY CURVE FITS WAS FXCFFDFD AT A TEMPFrATURF OF , 1PF12.5 01140 

2) 01150 

FPR2— 1 • 0 01160 

60 TO 355 Oil 70 

205 T5=I5+1 Cl 1 80 

60 TO 201 01190 

206 TF(TX?T(JiI)-TXK(I5-1,I) ) 210.220.215 01200 

210 T 5=15-1 C1210 

60 TO 200 01220 

215 XKB( J.I )=XK< IS-1, I > + ( (XK ( T5» I )-XK (15-1,1) ) / ( TXK ( 1 5. I ) -TXK ( 15-1 ,1) ) 01230 

1)*(TX2T(J,I) -TXK ( 15-t » I ) ) C1240 

60 TO 230 01250 

220 XKB(J,1)=XK(I5,I) 01260 

230 TF ( 16—1 > 203,203,231 C1270 

231 TF(I6-LCP> 232,232,203 01280 

232 TF(TX2T(J,I)-TCP(I6»T) > 234,245,233 01290 

233 T6SI6+1 01300 

60 TO 231 a 01310 

234 TF(TX2T(J»I)-TCP( 16-1 , I ) ) 235,245,240 01320 

235 T6=I6-1 01330 

60 TO 230 01340 

24 0 CPB ( J , I ) =0PX ( T6-1 » I )+< < CPX ( 16,1 )-CPX( T 6—1 , 1 ) ) / ( TOP ( 16 , 1 ) -TCP ( T6-1 , 01350 

IT) ) ) * ( TX2T ( J, T >-TCP ( T6-1 , T > » 01360 

60 TO 280 01370 

245 0PB( J, I >=CPX( T6, I ) 01380 

280 OONTINUF 01390 

15=2 01400 

16=2 01410 

300 OONTINUF C1420 

310 TF(DMP) 355,355,320 01430 

320 WRITE ( 6 , 330 ) 01440 

330 F0RMAT(/1X»32H PROPFRTIFS OF ABLATION MATFRIAL/) 01450 

WRITE (6,335 ) 01460 

335 FORMAT ( /5X , 5HYK ( I ) > 9x » 5HCP 1 1 > ,9X , 6HRH0 ( I ) / ) 01470 

WRITE (6,340 ) ( YK( I) »cP( I ) ,RHO(I) , I=1,NP) 01480 

340 F0RMAT(2X,1PF12.5,2X, 1PF1 2.5.2X, 1PE12.5) 01490 

wRlTE(6, 345) 01500 

345 FORMAT (//1X.32H PROPFRTIFS OF BACK-UP STRUCTURE/) 01510 

WRITE (6, 347 ) 01520 

347 FORMAT </5X»8hxKR< J, I) , 7X , AHCPP ( J, I ) , 7X , RhRHOBX ( I ) , 7X, 7HFMFB ( T ) , 8x , 01530 

17HEMBR ( I ) , 9X , 6HDXB ( I ) / ) 01540 

no 350 1=1 »NMR 01550 

Kt=NPM(I) 01560 

no 349 J= 1 ,KL 01570 

WRITE (6, 348) XKR ( J , I ) , CPB ( J, J ) , RmOBX ( T ) ,fMFB ( I ) ,FMBR( T) »DXB( I ) 01580 

348 F0RMAT(3X,1PF12.5,3X,1PE12.5,3X,1PE12.5,3X,1PF12.5,3X,1PF12.5,3X.1 C1590 

1PF12.5) 01600 

349 CONTINUE 01610 

350 OONTINUF 01620 

355 rFTURM 01630 

FND 01640 
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sibftc abl nonon 

THIS subroutine DETERMINfS thf mass flow ratf FROM thf nooio 

ABLATING NODES 00020 

SUBROUTINE ABt ATF 00030 

C 00040 

DIMENSION TITLE (1?) , HF ADNG ( 12 ) t X T ONT ( 1 2 , 1 2 ) »TKC(20) f XKC(?n) , 00 050 

lrPC(20) f TKV(20) , XKV ( 20 ) ,TCPV(20) ,CPV<20) ,TIME(30o) ,QCON(3O0) , 00 0*0 

2OPADC300 ) , VEL( 300) # XNPM ( 1 2 ) $ NKPR ( 1 2 ) ,NCPB< 12) ,TXk (20, 1 2> . XK < 20, 3 2) D0070 

3, TCP ( 20 #1?) ,CPX< 20 , 12) t RHOBX (12) , XBM ( 1 2 ) ,FMFB(12) , EMBR < 1 2 > , HXX ( 1 ? ) OOOAO 

4 , GAPX ( 12 ) 9 FTFST ( 12 ) ,BTEST(12) , TFmDI ( 20 0 ) , TX 1 ( 20 0 ) ,TXp(200) # 00090 

5TX2T(lOt 12) fTHLK200) ,TUL2<200) r HX ( SO > , Tw ( 50 ) f IR ( 50 ) , TR1 ( SO ) , POtOO 

6 T P2 ( SO ) # TUL ( 50 ) , I EM ( SO ) , T Y ( 200 ) , A ( 20 0 ) , B ( 200 ) 9 C ( 20 0 ) , P ( 200 ) , P011O 

7r< 50) 9 RHO ( 50 ) 9 CP (50 ) ,DXP (12) ,VKB( 10,12) ,CPB( 10# Ip) ,XMOG(50) 9 P0120 

8YK (50) 1 Ab (10*12) * BB (10*12) ,CR(10,12) , OB ( l 0 , 12 ) 9 SR (10,12) » 00130 

9«P1 (10, 12) ,RR2< 10, 12) f H( 12) ,5(50) ,NPM( 12) 00140 

PIMFNSION TTi.lL ( 50 ) 9 RHO Y 1 (SO) ,PHOy2(50) ,0RH0(50) ,TCPC(20> 00150 

C 00160 

COMMON TKC,XKC» TcPC 9 CPC ,TKV,XKV, TCPV 9 CPV , XNPM 9 RHoRX » XRM» F M BB * 00170 

lFMFB,NKPB,NCPP,TXK ,X k , TCP , CPX , NPM , GAPX , FTFST , PTEsT , T£MDT , TX 1 , 001 AO 

2tV 2,TX2T,TUL,TUL1,TUL2, TP,IRl,IR2r A,B,C,P,S,R,AB # RB,CP#DB,SB. 00190 

3rP 1 , RP2 , T Y , RHOY 1 9 RH0Y2 , XMPG , RHO , CP , YK , XKR , CPB , DXq , DT $ XLOST 9 0020 0 

4TABL,TCHAR,TRFC,RH0V,RH0C,FPL0W,fMV/, EMC,H300,NKC,NCPC,NKV, NCPV, 00210 

5NP»NMP,NPRS,NPF,TFST?, TFMPI ,TX0,tFNV,HEN\/,FFNV,QlOSS,TLIM,TINT 00220 

COMMON II, 12, T3, 14,15, If . PIN , TNT , OX , XMT , TL t VL , BL ,DMP , FRR1 , ERR?, 00230 

IF RR3 » FRR 4 , HV , VPT ,CHArk,CHARC,ABLk , APLC, XmOC,H 00240 

C 00250 

XMT=0,0 00260 

L T NT = TNT 00270 

kI=NP 002 AO 

TF(PMP) rt,8,3 00290 

3 wP I TE ( 6 , h ) 00300 

5 «fORMAT(//jX,?QHMASS FLOW FROM APlATING N ODES//) 00310 

8 DO 200 K KT=LINT,NP 00320 

TF(IRKKl)) 11 ,11,12 00330 

11 TF(TX1 (KX) .LF.TArD go to 9 00340 

12 TUL1 (KI )=AMAX1 (T(JL1(KT) ,TVl (KT ) ) 00350 

TP1(KT)=1 003*0 

GO TO 20 00370 

9 TF(TX1 (Kl )-TAPL) 10,10,20 00380 

10 RH0Y1 (KI )=RHOV 00390 

GO TO 50 00400 

20 TF(TULKKT)-TCHAR) 40,30,30 00410 

30 RHOYKKI )=RHOC 00420 

GO TO 50 00430 

40 RH0Y1 (KT )=RHOV+(RHOV-RHOf )* ( (TULl (KI )-TAPl ) / ( TAP| -TCH AR ) ) 00440 

SO TF(IR2(Kjt)> 52,52,54 004SO 

52 TF(TX2(KI) .LE.TABL) GO TO 56 00460 

54 TUL2 t K I ) -AM Ax 1 ( TUL2 ( K T ) , TX2 (K T ) ) 00470 

IRp(KT)=l 004 AO 

GO TO 70 nouqn 

56 TF(TX 2 (KI)-TaPL) 60,*n,7n 00500 

60 PH0Y2(KI )=RHOV 00510 

GO TO 95 00520 

70 T F ( TUL2 ( K T ) — TOHAR ) 9n, 80,80 00530 

AO RH0Y2(KI )=RHOC 00540 

GO TO 95 00550 

90 RH0Y2 ( K I )=RHOV+(RHOV^RHOC)*( (TUL2(K!)-TARL)/(TABl-TCHAR) ) 00560 
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95 0RH0<KI> = C<RH0Yl<KI>_RH0Y2(Kn)/nT>*DX 00570 

IF(KI-NP) 97,96,96 00580 

96 nRH0(KI)=nRH0(KI>/2.n 00590 

go to 98 00600 

97 IF(KI-INT) 96,96.98 00610 

98 IF ( DRHO ( KI ) ) 110.120,120 n0620 

110 ORHOIKI )=0.0 00630 

120 xMT=XMT+DRHO(Kl > 00640 

XMDG(KI)=XMT 00650 

IF(DMP) 190,190,150 00660 

150 WRITE < 6 , 160 ) XMDG < K I > , DRHO (KI) ,RhOY2 ( KI ) , RHOY1 < K I > 00670 

160 FORMAT ( 1X,5HXMDG=,1PE12.5»2X,5HDRH0=, iPFl 2 . 5. 2X ,6HRH0Y2=, 1PE12.5.2 00680 

1X,6HRH0Y1=,1PF12.5> 00690 

190 Kl=KI-l 00700 

200 CONTINUE 00710 

RFTURN 00720 

FND 00730 
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SIBFTC OXID F0000 

F 0010 

THIS SUBROUTINE CALCULATES THF HFATING RATE DUE TO COMBUSTION F005>0 

IT IS ASSUMED THAT OXYGEN AND CARBON REACT TO FORM CO ONLY. F0030 

F0040 

SUBROUTINE OXTDATtXMQO, OOXID) F0050 

C F 0060 

OOXID=XMDO*4000.0/36n(1.0 F0070 

OOXIDrO.O F0080 

RFTURN F0090 

FND F0100 
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THIS SUBROUTINE DETErMINfS THE FORWARD TTM£ STEP TEMPERATURES 
RY SOLVING THE TRI-DlAGONAL MATRIX 
SUBROUTINE SWUFTt ArR f C»DtTfN,nMP) 

DIMENSION A ( ?00 ) f B ( 200 ) , C ( 200 ) f D( 200 ) # T ( ?00 ) • CP ( 200 ) , DP ( POO ) 
CP( 1)=C(1)/B(1 ) 
nP(l)=D(l)/B(l ) 
no 100 I= 2 #N 

CP< I >=C< I >/<B< I )-A(I)*CP<I«-l) ) 

DP(I) = (D<I)-A(I)*nP(T-l) )/(B(T)-AlI)*rP(M) ) 

100 CONTINUE 
T(N)=DP(N) 

NM1=N-1 

no 2on j=i fNmi 
T = N-J 

T< I )=DP<I)-CP< I >*T( I+t ) 

200 CONTINUE 

TE(DMP) 300 t 300 9 250 
250 wR I TE ( 6 f 2B0 ) 

2B0 FORMAT ( / / 1 X f UXHCOEFF T C IENTS CALCULATED By SUBROUTINE SWUFT//) 
wR I TE ( 6 f 270 ) 

270 FORMAT ( 6X f 5HCP ( I ) f IOXf 5HnP ( I) f 1 0 X f 4HT ( I ) / ) 

WRITE (6f 275) ( CP < I ) r DP ( I ) f T ( I ) f 1=1 f N ) 

275 fORMAT(PXf1PE12.5f2X, IPEI P.S.PXf 1PE12.5) 

300 RFTURN 
FND 


FOOOO 

F0010 

F0020 

F0030 

F00UO 

F0050 

FOOftO 

F0070 

Fooen 

F0090 

F010O 

Fonn 
F0120 
FO I 30 
F0140 
FOisn 
F01B0 
FO 170 
FO 1 AO 
F01S0 
FOPOO 
FO? 1 0 
F02PO 
FOP30 
FOP40 
F0250 
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IbFTC REC 600 Oft 

GOOlft 

THIS SUBROUTINE DETERMINES THF FRONT FACF LOCATION AND CHAR MASS 600PO 

RFMOVAL RATE G0030 

G0040 

SUBROUTINE RFCESs ( XMDC f XI OST , TRFC f DT , RHOc * TS » SR t TX2 1 NREC . NRS # FRRs 9 G0050 

lsXO t SHOT , DMP ) GO 060 

C G0070 

DIMENSION TS(SO) »SR(sn> 6008ft 

TF < TX2-TREC ) 10,20*2n 60090 

50 XMDC=ft*0 GO 100 

XLOST=0.0 G0110 

SDOTrn.O GO I PO 

GO TO 60 GO 130 

?0 IF(NRS-1 )25*25r2l GOlUft 

21 IF ( NRS-NkEC^) ?2t 22*25 60150 

22 TF(TXP-TS(NRS) ) 32*4n*3ft G016D 

25 wP I TE ( 6 f 26 ) TX2 G0170 

26 FORMAT < 1H0 * 75H THE RANGF OF THE SURFACE RECESSION TA0I E WAS EXCFfD GOlflft 

1FD AT A TFMPERATURE OF *1PE52*5) G0190 

FPR5=1*0 G0200 

GO TO 60 G021 ft 

30 NRS=NRS+1 G0220 

GO TO 21 G0230 

32 TF(TX?-TS(NRS-1) ) 34 # 40 * 36 G0P40 

34 NRS=NRS-1 60250 

GO TO 2ft G0P60 

36 SX=SR(NRS-1 )+( (SR (NRS)-SR (NRS-1 ) ) / ( TS < NRs ) -TS < NRS-1 > ) ) G0270 

1*(TX2— TS( NRS— 1 ) ) G0280 

GO TO 5ft G0290 

40 SX=SR(NRS) G0300 

SO *1 OST=3ftO.O*SX*DT G0310 

XMDC= < XL OST ♦RHOC ) /DT G0320 

sDOT = SX*300 * 0 60330 

IF(DMP) 60 * 60 # 52 60340 

52 WR I TE ( 6 * 54 ) SX»XLOST # XMnc 60350 

54 FORMAT* 1H0*3HSX=, 1 PE 1 2.5, 3X#6HXLoST=tlPE 1 2.5* 3X,5HXMDC=M PEI 2*5) G0360 

60 RFTURN G0370 

FND 60380 
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IBFTC TEMP HOOOO 

THIS SUBROUTINE DETERMINFS THF INITIAL TFMPfRATUrE DISTRIBUTION HOOIO 

IN THE HEAT SHIELD STRUCTURE H0020 

SUBROUTINE TEMPO H0030 

C H0040 

DIMENSION TITLE ( 1 2 ) »HEADN6( 12 ) , XlDNT (12,12), TKC(20)»XKC(20> * H0050 

1CPC(20) »TKV(20) ,XKV(20> , TCPV ( 20 ) ,CPV ( 20 ) , TIME ( 300 ) * QCON ( 300 ) , H0060 

20PAD( 300 ) » VEL <300) . XnPM( 1 2 ) ,NKPB< 12) ,NCPB(12) , TXk ( 20 . 1 2 ) . XK ( 20, 12 ) H0070 

3. TCP (20 #12) ,CPX(20,12) ,RH0BX(12) ,XBM<12) ,EMFB(12) #EMBB(12) #HXX(12) HOOBO 
4,0APX(12) #FTEST(12) »BTEST(12> #TEMDI (200 ) . TX1 (200 ) ,TX2I200) . H0090 

5TX2T (10,12) #TllLl(200) ,TUL2(200) #HX(50) »Tw ( 50 ) » IR < 50 ) . IR1 ( 50 ) , H0100 

6lP?(50) ,TUL(50) »IEM(SO) ,TY(200) # A (200) #B(200) , C < 200 > ,D(200 ) , HOI 10 

7R ( 50 ) #RHO ( 50 ) ,CP(50) ,DXP(12) #XKB<10#12) #CP8(10#12) #XMD6(50) • H0120 

8YM50) »AB( 10,1 2) ,BB( 10,12 ) ,CR< 10.12) ,DB( 10,12 ), SB (10, 12) . H0130 

9RB1 (10,12) #RB?(10#12) , H( 1 2 ) ,S (50 ) , NPM( 12 ) H0140 

DIMENSION TTuL(50),RHOY1(S0) ,PH0Y2 ( 50 ) »DRHO ( SO ) ,TCPC ( 20 ) H0150 

C H0160 

COMMON TKC»XKC,TCPC, CPC. TKV.XKV, TCPV, CPV , XNPM , RHoBX , XRM , FMBB , H0170 

1FMFB , NKPB , NCPB , TXK , XK , TCP ,CPX ,NPM , GAPX ,FTEST , RTEsT , TEMDI . TX1 , HO 1B0 

2tX2,TX2T, TUL. TUL1 ,TUl2,IR»IR1»IR2»A,B»C,D»S»R»AB,BB,CR»DB*SB, H0190 

3RB1 ,RB2, T Y,RhOY1 ,RH0Y2, XMDG,RHO# CP , YK, XKR,CPB, DXr ,0T , XLOST • H0200 

4TABL ,TCHAR , TREC , RHOV . RHOC ,FBLOW , FMV , EMC ^H300»NKC, NCPC , NK v »NCPV» H0210 

5NP , NMR , NPBS , NPF , TEST? , TFMP I , TXO , TENV , HENV , FFNV , OLOSS , TL IM , T I NT H0220 

COMMON I 1,12,13, 14, IS, 16 ,OIN, INT , DX, XMT »TL , VL , BL.DMP, FRP 1 , ERR2 , H0230 

1FRR3 , FRR4 , HV , VPT , CHARK , CH ARC , ABLK , ABLC , XMDC , H H02U0 

C H0250 

x=0.0 H0260 

TF (TEST2 ) 300,100,200 H0270 

100 DO 150 L=1,NPF H02B0 

TX1(L)=TEMPI H0290 

TX2(L)=TEMPI H0300 

TUL1 (L)=TX1 (L) H0310 

TUL2(L)=IX2(L) H0320 

TEMDI (L ) r TEMP I H0330 

ISO CONTINUE H0340 

DO 160 I=),NMR H0350 

jNrNPM ( I ) H0360 

DO 155 M=1 ,JN H0370 

TX2T(M,I)=TEMPI H03B0 

155 CONTINUE H0390 

160 CONTINUE H0400 

GO TO 320 H0410 

200 DO 220 L=1,NP H0420 

TFMDI (L)=TXO+( (TENV-TXO) /TL)**X H0430 

TX1(L)=TEMDI(L) W0440 

TX2(L)=TX1 (L) H04S0 

TUL1(L)=TX1(L) H0460 

TUL2(L)=TX1(L) H0470 

X=X+DX H0480 

220 CONTINUE H0490 

L=NP+1 H0500 

no 270 1=1 ,NMB H051P 

KJ=NPM(I) H0520 

DO 250 J=1 ,KJ H0530 

TEMDI (L ) =TXO+ ( ( TENV-TXO ) /TL ) **X H0540 

TX1(L)=TEMDI(L) H0550 

TX2 ( L ) =T£MDI ( L ) H0S60 
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TV2TCJ,I)=TEMni(L) H0P70 

x=x+dxbu> mosso 

l.=L+l H0590 

2S0 CONTINUF HOfenn 

X=X+(GAPX<I>/12.0) H0610 

270 CONTINUE H0620 

00 TO 320 HOftSO 

AM ARRRITARY TEMPERATURE OISTPIPllTION CAN BF BEAn IN FROM INPUT W 06 U 0 

DATA IF TEST2 IS A NFGATTVF NUMrFR H0A50 

300 wPITE ( 6» 31 0 1 H06.B0 

310 FORMAT ( 1 HO » 79H THF VALUF OF TFSTj> WAS NFGATTVF, SUBPOUTINF TFMPO S M0A70 

IhOULD not havf REFN CALLFO.) HOABO 

FPR1=1.0 H0690 

320 rFTURN M0700 

pND H0710 
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IbFTC D0N2 toooo 

this supkoutine dftfrminfs thf tfmperaturf of points a fixed room 

distancf from a rfferfncf plane from the temperatures calculated ino?n 

in a varying thickness tooso 

tnnun 

e;l IBROI ITI NF Df)N2 ( XLOST t XARR AY . T ARRAY > NA . XnODE * TEMp , XNonEV , KK * XLSTV > IOOSO 

lnx) 10060 

C T0070 

DIMENSION XArrAY(SO) ,TAPRAY( 50) » XNODE « 50) »TEMP(5o) .XNODFVI50) I 00 AO 

C 10090 

K=0 T0100 

nXT=0.0 T011O 

do ion i=i »na 10120 

IF(XLSTV.LE.DXT) GO TO ISO T 0 1 30 

K=K+1 T0140 

100 OXT=DXT+UX T0150 

ISO kKsNA— K T0160 

XK=K 10170 

xNODEV ( 1 ) =XLSTV I01A0 

TFMPM >=TARRAY( 1 ) T 0 1 90 

DO 200 I=t»Kx 10200 

vNODE ( I ) =XK*DX-Xl_STV T 02 1 0 

CALL DISCT3(XN0DE( I ) , XARR AY »TARR AY , NA , TFMP < 1+1 ) ) T0220 

xNODEV ( 1+1 ) =XK*DX 10230 

200 XK=XK+1.0 T0240 

rFTURN T02SO 

FND T0260 
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SIBFTC UINTRP JOOflO 

9U8R0UTINF UlWTRP(X f yTBL f YrYTPLt N# J) JOOIO 

DTMFNFION XTBl <50) *YT*L<RO) J0020 

1 zj J0030 

T F ( I * GT . N . OR ♦T*LT*2) t = ? • J0040 

10 TF(XTPL( 1-1) .1 E.X.ANn.X.l F.XTPLCT) ) GO TO 40 J0050 

TF ( X . GT « XTBL ( T ) ) GO TO 30 J006O 

20 T - I — 1 J007O 

TF1I.GE.2) GO TO 10 JOOflO 

T = ? JO 090 

GO TO 40 J0100 

30 T”i +i jonn 

TF(I.LE*N) Go to 10 JO 120 

T=N J0130 

40 FPACT=(X-XTBL< I-l) )/(XTRL<I)-XTPL<I-l) ) J0140 

Y=YTBL(I-1 > + ( YTRL(I)-YTBi <I-1))*FRACT J0150 

PFTIJKN J01A0 

FND JO 170 
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SIBFTC TSOT K0000 

FUBROUTIt'JF ISOTHM ( DEPTH* TFMP , RONn» Na ANS ) K0010 

DIMENSION DEPTH ( 1 ) a T fMP( 1 ) KOOPO 

ans=-i. koo?o 

K-N— i Kooun 

no ion t=i»k koobo 

TF ( TEMP ( 1 ) —BOND ) ? 1 1 1 3 K0060 

1 ANS=OFPTH(I> K0070 

fiO TO 100 KOOflO 

2 jF<TEMP(I+l)-ROND)10n, 100.4 K 0090 

4 AMS=DFPTH <1 + 1 > - (TFMP ( X + 1 ) -BOND ) * ( DEPTH < 1 + 1 ) -DEPTH ( I ) ) / < TFMP < I +1 > - K 0 1 0 0 

lTFMP(T)) KOI in 

fiO TO 100 K0120 

3 TF (TEMP ( 1+1 ) -ROND )5»100»100 K0130 

5 ANS=(TEMP(I)-R0ND)*(DFPTH<X+1)-DfPTH<I) )/<TFMP(I)-TEMP(I+1> 1+DEPTH KOI 40 

id) KOisn 

100 OONTINUF K0160 

TF(B0ND.E(3.TfmP(N) ) An|9=DFPTH(M> K0170 

RFTURN KOI 80 

FND K0190 
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FUBROUTlNF SAVE ( S AVE l , SA VF2 » SA VF3 * USE , NXl , V ALUE » DT . TFTN AL » T I MF . 
1THIN6) 

niMFNSION SAVF1 ( 1 > ► SAVE? 1 1 ) *SAVF3< 1) 

I iSE— 0 • 0 

<;AVEl(NXl)=VAtUE 

NX2=NX1-1 

TF(Nx?.FQ.0)NX2=3 

SAVF2INX2) — VAI.UE 

NX3=NX2-1 

TF(NX3.EQ.0)NX3=3 

SAVF3(NX3)=VA| UF 

TF( (TTME.LT. <?.*DT) ) .OR. (TIMF.GF. ( TFINAL-3 .*DT > ) ) GO TO 4 
GO TO <1 »?»3) »NX1 

1 TF ( ( ( ABS<SAVF2«1 )— SA VF2 < ? ) n .LE. .001) .OR . < ARS < SAvE2 ( ? > -SAVE2 ( 3) ) 
l.t E. .001 ) )GO TO 5 

tF ( ( (SAVE211 ) ,LT.SAVf?(2> ).ANO. ( SA VF2 ( 2) .PT. SAVE? < 3) ) > .OR. ( (SAVF?< 
11 ) . T.SAvF 2<2) ) .AND. (SAVF?<2> ,LT,SAVE?(3) ) )) USE=sAVE 2 ( 2 ) 

5 THING=SAVF2(2) 

GO TO 4 

2 TF< ( (ABSISAVF3(1)-SAVF3(?) ) ) .1 E. . 00 1 1 . OR . t ARS t S AvE3 ( 2 ) -S AVE3 ( 3 ) ) 

1 .1 F. .001 ) )GO TO 6 

TFH ( SAVE 3(1) , LT , SAVF3 ( 2 ) ) . AND. ( S AVF3 ( 2 ) . GT . S AVF3 ( 3 ) ) > .OR. ( (SAVF3C 
11 > .GT.SAVF3I?) ) .AND. (SAVF312) ,LT.SAVE3(3) > > > USE=sAVF3 < 2 > 

6 THING=SAVF3(?) 

GO TO 4 

3 IF ( ( ( ABStSAVFI (1 I-SAVFK?) > ) ,I.E, .001) .OR. < ARS < SAvFl ( 2 > -SAVE1 (3) ) 

I .1 E. .001 ) )GO TO 6 

TF( ( (SAVE) (1) .LT.SAVF1 (2) ) . ANO. (SAVE1 (2) .GT.SAVF) (3) ) ) .OR. ( (SAVFl ( 

II ) .GT.SAVFK?) ) .AND. (SAVF1 (2) , LT . SAVE1 ( 3 ) ) ) ) USE=SAVF1 < 2 > 

7 THING=SAVF1 (2) 

4 MX1=NX1+1 
TF(NX1.EQ.4>NX1=1 
RFTUHN 

fnd 


L 0000 
l 0010 
10020 
L0030 
t 00UO 
10050 
10060 
L0070 
10080 
l 0090 
10100 
L 0 1 1 0 
L0120 
L 01 30 
L014O 
1.0)50 
I 0150 
l 0170 
L01RO 
10190 
L0200 
L0210 
10220 
10230 
10240 
10250 
L0260 
10270 
10280 
10290 
10300 
10310 
L 0320 
L0330 
10340 
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SIBFTC PISCT3 

FUBROIITTNF D I FCT3 ( X A , TABx »T ABY»Ny t ANS) 
DIMENSION TA8X(1) »TARY(1 ) 

CALL OISSFR < XA * TARX * \ »NY,?» NN ) 

NNNr3 

CALL LAGkAN(XA*TARX(nN) ,TABY<NN> ,NNN» ANS) 
RFTURN 

fnd 


moooo 

M0010 

M0020 

MQ030 

M004n 

M0050 

Mooen 

M0070 
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SIBFTC DISS WOOOO 

SUBROUTINE DISSER (XArTABt I rNXr ID»NPX) W0010 

DIMENSION TAB (2000 ) W0020 

C DIMENSION TAB (2000 > W0030 

NPT=ID+1 W0040 

NPB=NPT/2 W0050 

NPU=NPT-NPB W0060 

IF (NX-NPT) 10 r 5r 10 W0070 

5 NPX=I W0080 

RETURN W0090 

10 NLOW=I+NPB W0100 

NUPP=I+NX-(NPU+1) W0110 

DO 15 II=NLOW»NUPP W0120 

NLOC=II W0130 

IF (TAB(II)-XA) 15 r 20 » 20 W0140 

15 CONTINUE W0150 

NPX=NUPP“NPB+1 W0160 

RETURN W0170 

20 NL=NLOC-NPB W0180 

NU=NL+ID W0190 

DO 25 JJ=NL t NU W0200 

NDIS=JJ W0210 

IF ( TAB ( JJ ) -TAB ( JJ+1 ) ) 25>30»25 W0220 

25 CONTINUE W0230 

NPX=NL W0240 

RETURN W0250 

30 IF (TAB(NDIS)-XA) 40>35>35 W0260 

35 NPX=NDIS-ID W0270 

RETURN W0280 

40 NPX=NDIS+1 W0290 

RETURN W0300 

END W0310 
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SIBFTC LAGR 

SUBROUTINE LAGRAN (XA»X* Y *N> ANS) 
DIMENSION X (200 ) * Y (200 ) 

C DIMENSION X (200 ) * Y (200 ) 

SUM=0.0 
DO 3 1=1 t N 
PROD=Y ( I ) 

DO 2 J=1#N 
A=X(I)-X(J) 

IF (A) 1*2*1 

1 B=(XA-X(J) >/A 
PROD=PROD*B 

2 CONTINUE 

3 SUM=SUM+PROD 
ANS=SUM 
RETURN 

END 


T0000 

T0010 

T0020 

T0030 

T0040 

T0050 

T0060 

T0070 

T0080 

T0090 

T0100 

T0110 

T0120 

T0130 

T0140 

T0150 

T0160 



SIBFTC MORE NOOQO 

DIMENSION TITLE* 12 > »X(2000) »Yl (2n00 ) . Y2 (2000) , Y3 ( 2000 ) .Y4(2000) N0010 

REWIND 11 N0020 

RFAD(ll) (TITLE(I) *1=1*12) N0030 

RFAD(11)X(1) ,Yl(l> »Y2(1> * Y3 ( 1 ) * Yto ( 1 ) N0040 

Y3(1)=Y3(1)*12.+Yl(l) N0050 

Y4(1)=Y4(1)*1?.+Y1(1) NOOAfl 

1=2 N0070 

30 RFAD(11)X(I)»Y1(I)»Y?(I),y3(I)»Y4(I) NOOSO 

TF(X(T)-5001.)10,20,zO N0090 

10 Y3(I)=Y3(I)*12.+Y1(I) N0100 

Y4(I)=Y4(I)*12.+YMT> N0110 

T = I+1 N0120 

GO TO 30 NO 1 30 

20 NPL0T=I-1 NO 140 

YMISYKI) N0150 

YM2=Y2(1> N0160 

YM3=Y3(1) N0170 

YM4=Y4 ( 1 ) NOlflO 

DO 40 K = 2 .NPLOT N0190 

IF <Y1(K).6T.YM1> YMi = Y1(K) N0200 

IF (Y?(K> .GT.YM2) YM? = Y?(K) N0210 

IF <Y3(K) .GT.YM3) YM3 = Y3(K) N0220 

IF (Y4(K) .GT.YM4) YMn = Y4(K) N0230 

40 CONTINUE N0240 

1000 FORMATUHl , (1?A6) ) N0250 

CALL ACC END (X.Y1.Y2.Y3.Y4 .NPLOT > N0260 

XMAX=X (NPLOT ) N0270 

CALL APLOT (X.Yl .XM aX.YMI .TITLF, NPLOT) N02RO 

CALL RPLOT (X,Y2,XMAx, YM?, TITLE) N0290 

CALL CPLOT (X.Y3.Y4.XMAX, YM3 * YM4 « T ITLF • Y1 ) N0300 

WRITE (6* 1000) (TITLE! I) .1=1*1?) N0310 

WRITE (6. 1001) (X( I ) .Yl ( I) ,Y2(I ) ,Y3( I) »Y4( I ) . 1 = 1 * NpLOT) N0320 

1001 FORMAT (5E20. A) N0330 

wR I TE ( 6 , 1 0 02 ) XM AX » YM 1 , YM? , YM3 , YM4 . NPLOT N0340 

1002 FORMAT ( ///6H XMAX=F1 0 . 4, RH YM1=F10.4,RH YM2=F10.4»5H YM3=F10.4.5h N0350 

1YM4=F 10.4 .2XAHNPL0T=I4 ) N03AO 

RF AD (11) (TITLE (I),t = 1.12) N0370 

RFADdl)X(l) .Yl(l) »Y2(1) .Y3(1).Y4(1) N03A0 

1=2 N0390 

IF(X(1)-5001. )30, 50,50 N0400 

50 WRITE(6»1003) (TITLE(x) ,1 = 1 .12) N0410 

1003 FORMAT (////12A6) N0420 

RFTURN N0430 

FND N0440 
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SIBFTC ACCEN POOOO 

SUBROUTINE ACCEND (X»V*A»B»C»N) POOIO 

DIMENSION XlllfYllliAllI.RIDiCdl P0020 

K=i P0030 

101 SMALL=X(K) POOUO 

DO 100 I=KrN P0050 

DUMYsX(I) P0060 

sMALL=AMiN1(SMALL»OUmY) P0070 

if(Small.eq.x(I))INDex=i poobo 

100 CONTINUE P0090 

X ( INDEX ) sx (K ) P0100 

X ( K JsSMALL P0110 

SAVEsY(K) P0120 

Y(K)SY(INDEX) P0130 

Y(INDFX>sSAVE POltoO 

SAVEAsAOO P0150 

A (K )SA ( INDEX ) P0160 

A ( INDFX ) sSAVEA P0170 

SAVEBsB(K) POIBO 

B ( K ) SB ( INDEX ) P0190 

R( INDEX )sSAVEB P0200 

SAVECsC(K) P0210 

C (K )-C ( INDEX ) P0220 

C ( INDEX ) eSAVEC P0230 

KSK+1 P0240 

IF(K,FQ.N) RETURN P02S0 

BO TO 101 P0260 

FND P0270 
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*IbFTC AP LOT 

CJIHROUTINF AP» OT ( X, Y, XL TM,YLTM,TlTt E, IP| OT) 

0 T MFN4 I ON X< 300) , YTTTLE( 1 0) ,XTITLF( 1 0) 
nT MENTION TITL E( IP) » Y( 300) » ALoNGy (7) 

COMMON /mPC / AtL0W(7) , A| ONGX<7) , NPLOT , 7FR0 , X* AX , TF I V 

f ) A T A (XTITLF(T) , I=1,10)/34H TIME < 9EC . ) / 

HA1 A <YTiTLF(T).I=l»in)/3AH SURFACE RFCPS^IOn (IN.) / 

7FRO=n.O 

A I LOw(l)=^0. 

Al L Ow ( 2 ) r 1 00 • 

A I LOw(3)=PbO. 

Al i.Ow(4)=FOn. 

Al L0w(5)=l 00n. 

Al LOw<6)=PbOO. 

Al LOw ( 7>=F00o. 

MPL.0T = IPLOT 
HO 10 1=1,7 

1 T = T 

TF(XLTM- ALLOW ( I ) ) P 1 1 1 2 0 » 1 0 
10 PONT I NIJF 

30 WRITE (F,1000> XlTM,y1IM 

1000 FORMAT ( / // 77H APL.OT CANNOT PF DOnF PECAHFF FTTHFr XL T m F XCEFOFD FO 
Ion. OR YLTN h yCFf OED F. /FH XI IMrF IP . F , Fx , 6H YL T M=E Ip.*) // 1 9H WF 
?i\)OW OO TO BPLOT /// ) 

RETURN 

PO yMA X = ALLOW (it) 

|RIX=TI 
no ao i = .i ,4 
tT = I 

1 F ( YL T M *100. -AlLOW(T) )F0,F0,4n 
UO fONTinuF 
oO TO 30 

FO Y^AX “ALLOw (IT) /100. 

CALL RSTFRM 

CALI OKTUCN ( 123.1 OP ^,?4,Q24, 18,1 4,5, F) 

CALL PLOI1 d , 1 ,ZFKO.VMAx,ZFRO,YvAX,X,Y,nPLOT,1 , 1 W/ > 

Al ONGY( 1 )=0.U 
Al ONGY( 1 )=0.O 
no GO 1 = 1 , b 

CALI- LARtl X ( ALONPX ( T 1 , 1 ) 

CALL LABEl Y < ALONGY < I ) , 1 > 

AlOMGY(T+1)= ALONOX(T) + .P* XVAX 
FO A I ONG Y (1 + 1 ) = ALONGY(t) + .o* YM AX 
r A| L PRTNT(?OO,P7S,1p,0,3«,XTTTl .f ) 

CALL PR T .nT ( 4 7 , 200 , 0 , 1 ? , 3rt , YT T TLF ) 
f'ALL PRT-yT ( IPX, 1 00 0, 1 P, 0p72,TTTI f ) 

PALL OMPbliE 

wFTUKN 

FND 


onnno 

00010 

oonpn 

00030 

onnun 

oonsn 

ooofo 

O0070 

00040 

onngn 
ooi on 
oni in 

ooi po 

ooi 30 
oni un 
ooi bn 

OOI Fn 
00170 
00140 
00190 

onpno 

OOP 1 0 
onppo 

O0P3O 

O0P40 

oOPfSn 

OOPFO 

O0P7O 

O0P4O 

onpqo 
00300 
00310 
003PO 
00330 
00340 
003FO 
003FO 
00370 
00340 
00390 
00400 
OHM 0 
004 PO 
00430 
00440 
o04*sn 
ooufo 
00470 
00440 
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SIBF7C PPL OT ROOnn 

PiltpPOMTT NF PPI OT (X» y» XL T M r YL TM p T T TLE ) POOIO 

dimension X(3D0) 9 y ( 3nn) » yt itlf (in) , alonoy(7> , yt itl E ( l n ) Ronpn 

n T MFN9 I ON TlTlE(lP) R0030 

rOMMQM /aPC / ALLOW ( 7 ) p Al ONGX ( 7 ) , NPL OT pZFROpXMAX, TFI X ROOUO 

H At A (XT1TLF { T ) f 1 = 1 f 1 (7 ) /3PH TIME (SFC.1 PDOBO 

nAT A (YTlTLF(T)f I = lfin)/3«H PONDLINF TFmPERATIJRF (R> P 00BO 

Al 0NGV( 1 ) =0 , 0 R0070 

no to i = i # 7 Ronflo 

tT=i ponqn 

tF(yltm —allow { i ) ) 20*20,10 poioo 

io roNTiMiJF ROttn 

IVPITE (5*1000) YLTM P0120 

1U00 f-ORMAT(/// 37W PHLOT WILL NOT BF DONE BFCAUSF YLtMz F12.5 ////) P013n 

pftukn pomn 

P0 YPAX =ALlOW(IT) POIfSO 

tall rstprm ROifsn 

call ORll'ON ( l?3f 1 023, ?4, P?4» 1 tA* 5,5 ) P0170 

CALL PLOfl ( 1 » 1 p ?FRO , XMAY ,ZFRO,YmAX,X* Y 9 NPLOTpI, 1H/ ) R01PO 

HO 30 1 = 1,6 PO 190 

CALL LARH X ( ALONGX ( r ) , 1 ) P0P0O 

rALL LAPF:I Y ( ALONGY ( T ) * 1 ) ROPIO 

30 A| ONG Y ( I + 1 ) = ALONGY(T) + .2* YM AX R0P2O 

CALL PRINT (20 0*975* 1 Pp Op 3 «pXTTTLF) P0230 

CALL PRInT(47p20Dp0p IPp3ApYTITLF) RQPun 

CALL PRINT ( 123, lunO r 1 Pp 0 ,72 pTTTL F) P0P50 

CALL nMPhUF POPftO 

RFTLIKN R0P7O 

f- NO P0280 
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sibFTC CPLOT qonon 

FI iRROt IT T imF CPI OT ( X , Y 1 , Y? , XL IM # Y| TMl , YLI M2 * T I TLF , Y) qnflin 

dimension X(^nO) • Yl ( ^no) • Y2(300 ) . ytiti E( i o) fYY(?ono) *\TTtle:( in) foopo 

mimfnfio/m Tin Ed?) #v(3nn) >alongy(7) fooso 

MTMFN^ION CURVE (1) fVpllG<u>»HP1'G(7) F0040 

COMMON! / ARC / ALLOW < 7 ) > A| ONGX ( 7 ) , MPI.OT r 7FP0 t XMAX „ TFI X F00SO 

nATA (VRuo(j) , I = d4) / mn.o»fin.n» 20 .n»)n.o / fooso 

OAT A (HRUG(I) , 1=1*7) / 1 . 0 * 2 . n , fS t 0 $1 0 . 0 , ?0 . 0 , SO . n *1 Of) . 0 / F0070 

hAlA (XT 1TLF ( T ) * 1 = 1 r 1 0 ) / .38H TIME ( FFO . ) / FOOpn 

HATA ( YT jl TLF (T)*I= lrt0)/3OH nTSTAMCF (IN.) / qonqn 

PAT A ONF/uHIofO / * T Wo/ 4M 1 46 0 / FOtOO 

data won/ i hi /»too/ih? / foiio 

C *** FOUR (4) CHARACTERS ARE ALLOWFD FOR CURVF(l) FQ1?0 

Cl'RVE(l)=ONF S013n 

HFAC7 RrHrtHG ( T F I X ) F0140 

FYMPOL=WuN FOlsn 

YFIG 3 AMaX 1 (vLlMnYl.TM? ) FOlbO 

NiCURVF =1 F017P 

no 1 T=1fMPL0T FOIflO 

1 Y Y ( I ) “ Yi(I) FOlQO 

no 7 t=i »u fo?oo 

tt= I Fo?in 

T F ( Yb T G+ 1 DO • - Al LOW ( T ) ) 6 # 6 , 7 FO??0 

7 CONTINUE F0P.3O 

WRITE (6HUflll) YtTMl.YLTN? F0?40 

1000 FORMAT (/// 3QH (.PLOT WTH NOT PF DONE PFCAUFF Y[ IM 1 =F1? . F , ) pH OR Fppsn 
1 V I lM2r F 1 ? . 5 //// ) fo?60 

RF TURN F0P7P 

6 YMAX -ALLOW (Tl)/100. FP28P 

VFAC I R-Vb( 'G ( I T ) FO ?90 

CALL RSThPM F03n0 

CALL GRTtjr-N ( 123# 1 0? } * 24 . Q24 » 18# 1 8 * S » b ) 

0=1 F03pP 

70 no 10 I -J * NPloT fo33o 

TI= I F034P 

TN I 3 T — 1 C038P 

|F( YY(n-Yd) >20*1/1*10 C0XAO 

10 C0NT1MUF F037P 

MOPT3NPLOT-J+1 F03AP 

t I 3 J + NOPT/2 Fp3QP 

TVLOC3(YmAX-YV(LL) >* 18 .* w FACTP -»-?4 * -4. F0400 

THL0C3 X(LL)*l8. /HFncTP +123. -48. cpam 

CALL PRTi\jT(IhI 00, TVLOC* 8 * 0 * 4 , Cl IR VE ) Fp4?n 

CALI, p L0ri (1 , 1 *7FR0* XMAX,7F.RO r YMAXrX( J) , YY( J) , NOpT # 1 f SYM^Ol ) F043O 

TF (NCUPVF-1 1 Pnr8E f Q0 F 04 4 n 

8b no 88 I=1»NP| nT CQUEn 

8fe YY(I)3Y?d) FOUbn 

CHRVEd >=TWO F047P 

bYMPOL=TOO Fp4fln 

KjCURVF = p FpUQP 

J31 fpspp 

GO TO 70 FO^in 

20 NPT 3 I T — J F0S?P 

I l = J + i\PT/p FOS30 

IVL0C3(YMAX-YY(LL) )*18.^\/FACTR +24. -4, F0F4P 

T PL OC" X ( I L ) * 1 8 » /HFaCTP +123. -48. FOSSP 

CALL PRTNT(IhI OC»TVLoC* A, 0*4, CURVE) FQSbP 
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CALL PLOM (1,1 *?FRO*xWAy ,7ERO*YMAX t X( J) ,YY( J> , NPT . 1 • 9YMP0L ) 90570 

no 50 IJ= II* NPLOT 505*0 

JJz IJ 50590 

TP ( YY ( I J ) — Y(Tj) )50 # uo, 40 505no 

50 CONTINUE C051O 

IF(NCHRVt-l) Q0*ft5*P»l 50520 

40 jz jj 50550 

GO TO 70 50540 

90 Al ONGY( 1 )rO*n 50550 

no ion T = i,5 c 065n 

CALL LABFl X(AI ONGX(I) * 1) 90570 

CALL LARtt Y C Al ONgY(T) #1) 505*0 

100 Al ONbY(I + l )=Al ONGY(T) + ,?*YMAX 90590 

CALL PRINT (20 0,975* 1 2 , 0 , 5 A , XT T TLF ) 9070 0 

CALL PR I NT ( 47 , 20 0 , 0 , 1 2 , 3ft , YT I TLF ) 90710 

CALL PRTNT( 125, 1 000 , l?, 0,72,TTTLF) 907?O 

CALL nMPbllF 90750 

RFTUKN 90740 

f 5'D 90750 



APPENDIX D 


FORTRAN 

A 

AB 

ABLC 

ABLK 

B 

BB 

BL 

BLTEM 

BTEST 


C 

CB 

CHARC 

CHARK 

CP 

CPB 

CPC 

CPV 

CPX 

D 

DB 

DELTT 


PROGRAM TERMINOLOGY 

Description 

"A" coefficient in matrix, single subscript 
"A" coefficient in matrix, double subscript 
specific heat of material at TABL 
thermal conductivity of material at TABL 
"B" coefficient in matrix, single subscript 
"B m coefficient in matrix, double subscript 
Total thickness of backup structure 

value of 1460 isotherm depth from previous time step 

test to determine mode of heat transfer out of back surface of 
backup materials 

"C" coefficient in matrix, single subscript 

"C 11 coefficient in matrix, double subscript 

specific heat of material at TCHAR 

thermal conductivity of material at TCHAR 

specific heat of a node in ablation material 

specific heat of backup material node 

specific heat values in char specific heat table 

specific heat values in virgin specific heat table 

specific heat values in backup material specific heat tables 

"D" coefficient in matrix, single subscript 

"D ,T coefficient in matrix, double subscript 

time step in the time step table 



FORTRAN 


Description 


DMP 

DRH0 

DT 

DTS 

DX 

DXB 

DXV 

DXX 

EMBB 

EMC 

EMFB 

EMV 

EMX 

END 

ERRl ^ 

ERR2 

V 

ERR3 

ERR4 j 

FBL0W 

FC0NV 

FENV 

FRAD 

FTEST 


test used for dumping (DMP = 0 skip dump, DMP = 1.0 start dumping) 

local mass flow rate of ablation gas 

time step from the time step table in hours 

time step from time step table in seconds 

thickness of a node in the ablation material 

thickness of a node in a backup structure material 

variable ablation node thickness 

fixed ablation material node thi 

emissivity of back surface of each material in backup 
char material emissivity 

emissivity of front surface of each material in backup 
virgin material emissivity 

emissivity of front surface of ablation material 
code word for plot routine 

Control numbers for printing error statements when an input or 
calculational mistake is made 



blowing efficiency in reducing convective heating 

factor to correct convective heating rate for various body locations 
emissivity - view factor product to cabin interior 

factor to correct radiative heating rate for various body locations 

test to determine mode of heat transfer into front surface of 
backup materials 
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FORTRAN 


Description 


FV 

G 

GAPX 

H 

H300 

HEAD 

HEADNG 

HENV 

HTX 

HV 

HW 

HX 

IEM 

I PRC 

IPRCT 

IR 

IR1 

IR2 

NCASE 

NCPB 

NCPC 

NCPV 

NKC 


view factor for external environment 
defined by FORTRAN statement 
gap width between backup materials 
film coefficient between backup materials 
enthalpy of air at 300° K 

any 72 alphanumeric characters used to identify problems being 
run - printed at top of first page of output 

any 72 alphanumeric characters used to identify each input section 

film coefficient to cabin environment 

total enthalpy 

heat of degradation of virgin material 

wall enthalpy computed from enthalpy — temperature table 
enthalpy values in enthalpy table 

test used to determine if front surface is virgin or char for using 
proper emissivity 

variable print frequency in time-step table 
present print control number 

test to determine if node temperature is greater than TABL 
test used in determining node density at TX1 temperature 
test used in determining node density at TX2 temperature 
number of problems to be run 

number of points in each backup material specific heat table 
number of points in char specific heat temperature table 
number of points in virgin specific heat temperature table 
number of points in char thermal conductivity - temperature table 
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FORTRAN 


Description 


NKPB 

NKV 

NMB 

NP 

NPBS 

NPF 

NPL0T 

NPM 

NHP 

NPTT 

NREC 

NTRAPT 

NXA 

NXB | 

NXC ) 

NXD 

NXE 

J 

QBL0CK 

QC0N 

QC0NX 

QHW 

QIN 

QL0SS 

Q0XID 


number of points in each backup material thermal conductivity table 

number of points in virgin thermal conductivity temperature table 

number of materials in backup structure 

number of node points in ablation material 

total number of node points in backup structure 

total number of points in heat shield structure (NP + NPBS) 

output plot control number 

number of nodes per material in backup 

number of points in enthalpy — temperature table 

number of points in time -step table 

number of points in surface recession — temperature or time table 
number of points in trajectory input table 


dummy indices for subroutine SAVE 


amount of convective heat blocked due to mass injection into 
boundary layer 

trajectory table convective heating rates 
cold -wall convective heat rate at present time step 
hot wall convective heat rate without blowing 
net heat flux into front surface 

boundary condition for heat transfer to cabin interior 
heating rate due to combustion 
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FORTRAN 


Description 


QRAD 

QRADX 

QUIT 

R 

RBI 

RB2 

RH0 

RH0BX 

RH0C 

RH0V 

RH0Y1 

RH0Y2 

S 

SD0T 

SAVEIT 

SAVEXX 

SAVX 

SAVY1 

SAVY2 

SAVY3 

SAVY4 

SAVY1X 

SAVY2X 


trajectory table radiative heating rates 
radiative heat flux at present time step 
code word for plot routine 

thermal resistance due to conductivity between nodes in the ablation 
material 

thermal resistance due to conductivity between past and present node 
in backup material 

thermal resistance due to conductivity between present and forward 
node in backup material 

density of an ablation material node 

density of individual materials in backup 

mature char material density 

virgin ablation material density 

density of node at past time step 

density of node at present time step 

thermal capacity of a node in the ablation material 

surface recession rate 

depth of 1060 isotherm at any given time 
time corresponding to maximum depth of l460 isotherm 
time corresponding to maximum depth of 1060 isotherm 
surface recession depth at maximum 1060 isotherm depth 
bondline temperature at maximum 1060 isotherm depth 
term that will contain maximum depth of 1060 isotherm 
depth of 1460 isotherm at maximum 1060 isotherm depth 
surface recession depth at maximum l460 isotherm depth 
bondline temperature at maximum 1460 isotherm depth 
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FORTRAN 

SAVY3X 

SAVY4X 

SR 

T 

TABL 

TCHAR 

TCP 

TCPC 

TCPV 

TEMDI 

TEMPI 

TENV 

TEST2 

TDMP 

TIME 

TINT 

TITLE 

TKC 

TKV 

TL 

TLIM 

TKEC 

TS 

TTABLE 

88 


Description 

depth of 1060 isotherm at maximum 1460 isotherm depth 
term that will contain maximum depth of l460 isotherm 
surface recession values in surface recession table 
present time 

temperature at which ablation starts 
temperature at which ablation stops 

temperature values in backup material specific heat tables 
temperature values in char specific heat table 
temperature values in virgin specific heat table 
arbitrary initial temperature distribution values 
constant initial temperature distribution value 
interior cabin temperature 

test to determine proper heat shield initial temperature distribution 

time to start dumping or printing information used in checkout of 
program (sets DMP = 1.0) 

trajectory table time values 

starting time of problem 

control card used for reading in new data for successive problems 
temperature values in char thermal conductivity table 
temperature values in virgin thermal conductivity table 
total thickness of heat shield structure (VL + BL) 
time limit of problem 

surface temperature or time at which char removal is to start 
temperature or time values in surface recession table 
time values in time -step table 


J 



FORTRAN 


Description 


TTUL equals TUL if VPT = 0 or equals TX2 if VPT = 1 - used in computing 
char properties 

TUL maximum value of TX1 and TX2 

TUL1 maximum TX1 values - used in computing gas ablation rate 

TUL2 ■ maximum TX2 values - used in computing gas ablation rate 

TV sink temperature of external environment 

TW temperature values in enthalpy table 

TX1 temperature of nodes at past time step 

TX2 temperature of nodes at present time step 

TX2C temperature at fixed locations in ablation material as defined by XC 

TX2T temporary storage of TX2 temperatures for computing thermal properties 

1XK temperature values in backup material thermal conductivity tables 

TX0 initial temperature at front surface of heat shield for computing 

linear temperature gradient 

TY temperature distribution at forward time step 

VEL trajectory table velocity values 

VELX trajectory velocity at present time step 
VL initial virgin material thickness 

VLI initial ablation material thickness 

VLTEM value of 1060 isotherm depth from previous time step 
VLV variable ablation material thickness 

VPT test to determine if properties are irreversible with temperature 

WEKEEP depth of l460 isotherm at any time 

XBM thickness of individual materials in backup 

XC fixed location of nodes in the ablation material 


XI 


node number 



FORTRAN 


Description 


XIDNT 

XK 

XKB 

XKC 

XKV 

XL0ST 

XLSTI 

XLSTV 

XMDC 

XMDG 

XMD0 

XMDT 

XNP 

XNPM 

XPL0T 

XV 

YK 

YPL0T1 

YPL0T2 

YPL0T3 

YPL0T4 

ZZZ 


any 72 alphanumeric characters to identify each material 

thermal conductivity values in backup material thermal conductivity 
table 

thermal conductivity of backup material node 

thermal conductivity in char thermal conductivity table 

thermal conductivity value in virgin thermal conductivity table 

amount of solid ablation material lost in a time step due to surface 
movement 

distance from original surface to present front surface location, 
inches 

distance from original surface to present front surface location, 
feet 

mass loss rate of char 

mass gas ablation rate due to pyrolysis of virgin material 
mass flux rate of oxygen to surface 
total ablation rate 

number of nodes in ablation material 
number of nodes per backup material 
time to be written on tape and plotted 

location of nodes in variable ablation material thickness 

thermal conductivity of a node in ablation material 

recession depth to be written on tape and plotted 

bondline temperature to be written on tape and plotted 

1060 isotherm depth to be written on tape and plotted 

1460 isotherm depth to be written on tape and plotted 

ratio to determine when the limiting value of heat blockage has been 
reached 
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APPENDIX E 


GENERAL PLOW CHART 


Start 



Rewii 
Read I 

id 11 
TOASE 


Initalize 
program constants 


Read and vrite 
input data and 
setup initial conditions 


Transfer J;o subroutine 
TEMFD 

determines initial 
temperature distribution 


Calculate surface 
heating conditions 


Transfer to subroutine 
RECESS 


calculates surface recession 
depth and char ablation rate 


Transfer to subroutine 
C0EFF 

calculates tri-diagonal 
matrix coefficients 


Transfer to subroutine 
SWUFT 

solves tri-diagonal matrix 
for temperature distribution 


I Continue 


Transfer to subroutine 
PR0P 

calculates thermal 
properties 











92 





































TABLE I. - SAMPLE PROBLEM INPUT 
(t>) Fortran data card listing 


1 

TrPICAl ChARRJLMG ABLATOR - TEST Ca 6 E - 4/6/66 DONALD M. CURRY 


+ fo0U • 0 

+ 1)0 

+ 0 .0 


? 

1 







■MI*U 

+ u«« 

+ 1) . 1 

+on 

i on 








+ fcj 0 0 • u 

+00 

+ 0,1 

+un 

1 00 








+ 1.0 

+0 0 

+ 1.0 

+ 00 









TYPICAL 

CHmRHXnO 

A HI 

ATION matfrtai. propfrtifs 




+lMfoU . n 

+u 0 

+ i^bu . 0 

+ 1)0 

+ u • 0 

+ 00 

+ 34,0 

+ 00 

+ 20.n 

+ 00 

+ 0,00 

+ 00 

+ 0-65 

+ u 0 

+ 0 .76 

• +00 

+ 129.06 

+ 00 

+ 1.50 

+on 

+260.0 

+ 00 

+ 0.0 

+nn 

+ 1 . u 

+ u 0 

*M’,0 

+ U0 

+ 0.12 

+00 

+ 0.43 

+00 

+0.070 

+ 00 

+ 0.43 

+on 

31 

a 

2 Q 

2 ? 








+iubo . n 

+ U0 

+ 0.1? 

+ 1)0 

+ 1.0 

+ 04 

+ 0.1? 

+ 00 





+1460,0 

+ u 0 

+ 0.43 

+ u 0 

+ 1.0 

+ 04 

+ 0.43 

+nn 





+ 360 . n 

+ U0 

+ t» . OFb 

+ un 

+460.0 

+un 

+0.065 

+00 

+560.0 

+00 

+0.0656 

+ 0 n 

+ b6 0 * U 

+ 1)0 

+U .06b 

+ un 

+ 76 0 , 0 

+ 0 0 

+0.067? 

+00 

+860,0 

+ 00 

+0.0684 

+00 

+ 960 • u 

+ u 0 

+0.069 

+ 0 0 

+1060.0 

+un 

+0.U70 

+ on 

+ H60.0 

+ 00 

+0.070 

+ 00 

+360, u0 

+ 1)0 

+ 0.43 

+ 1)0 

+1100,0 

+on 

+ 0.43 

+ 00 





+ 0.U 

+ 0 0 

+ 9.0 

-04 

+bO0. 0 

+uo 

+ 9.0 

-04 





NO T K A J^CT OR Y - 1 

3=96 

PTH/Str-SOFT 






+ U . 0 

+ 1)0 

+ 96. n 

+ un 

+ 0.0 

+00 

+2.9?5 

+ 04 





+bPU • U 

+ U 0 

+w^.n 

+ u 0 

+ 0.0 

+00 

+2.925 

+ 04 





1 

3 

+ 0,1 

+ 00 









+ 3.0 

+ 1)0 











9 

2 











BACKUP ^ATbPIAL 



0 .] 

L INCHED 

THTCK 




+3P0.U 

+ 00 

+ 0 . 1)66 

+ 1(0 

+460.0 

+ 00 

+0.065 

+ 00 

+560.0 

+ 00 

+0.0656 

+ 00 

+ hNJ . U 

+ 00 

+0 . 06b 

+ 00 

+ 7 6 u • 0 

+ 0 0 

+0.U67? 

+ 00 

+860.0 

+ 00 

+0.0684 

+ 00 

+ 960 . 11 

+ 1)0 

+0,069 

+ 00 

+1060 . 0 

+ 00 

+ 0,07 

+ 00 

+1160.0 

+ 00 

+ 0.07 

+ 00 

+ 36U . 00 

+ 1)0 

+ 0.43 

+ U0 

+1100.0 

+00 

+ 0.43 

+ 0 0 





+ 34.0 

+ u 0 

+ 0.1 

+ 1)0 

+ 0.9 

+on 

+ 0.9 

+ 00 





+ O.0 

+ U" 

. .. . u 

+ 1)0 

+ U.U 

+00 

+ 0.0 

+ 00 





AT TKAN^FFh TO C AP I N 

F<\/VIROMMf- N f 

- HEMVrO.O 





+ 560 . (I 

+ 00 

+ 0 .0 

+ 1)0 

+ 0.0 

+ U 0 

+ 0.0 

+on 





IMTT1AL 

. |FMPhRATl)PL 16 

CONST AMT 







+ 0.0 

+ 1)0 

+630 . U 

+ 1)0 

+630.0 

+u n 
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+ u • 0 

+ nn 

+ 0 - n 

+ 00 

+342. 9 

+nn 

+1400.0 

+ 00 

+449.7 

+ 00 

+1800. 0 

+ 00 

+ t> 1 7 • 2 

+ u 0 

+2400. n 

+ un 

+791.0 

+ 00 

+3000 . n 

+ on 

+978.0 

+ 0 0 

+3600 . 0 

+ 00 

+ 1 1 1 3 , 0 

+ 11 O 

+4ouo . n 

+ 1)0 

+1200. 0 

+ 00 

+4224 . 0 

+ nn 

+1 300.0 

+ 00 

+4486.0 

+00 

+ 1 A 0 U . 0 

+ 00 

+ u n 

+ u0 

+1600.0 

+ U0 

+4936.0 

+on 

+1600.0 

+ 00 

+5127.0 

+ 00 

+ 1700 .0 

+00 

+6299.0 

+ 00 

+1600.0 

+ 0 0 

+5454.0 

+on 

+l°on.n 

+ on 

+5696.0 

+ 00 

t 2 . 0 0 0 • 0 

+ 00 

+6726.0 

+ 00 

+2100.0 

+un 

+6051.0 

+ 00 

+2200.0 

+ 00 

+5Q68.0 

+ 00 

+2^00.0 

+ 00 

+ b07ft. 0 

+ u 0 

+2400.0 

+nn 

+6186.0 

+un 

+2600.0 

+ 00 

+6291 .0 

+ 00 

-, 2600.0 

+ u 0 

+6395.0 

+ 00 

+ 2700. n 

+ 00 

+6497.0 

+00 

+2800.0 

+on 

+6697.0 

+00 

+^ Q i)0,n 

,+ QO 

+ t>699. 0 

+ 00 

+ 3000. 0 

+00 

+6005. n 

+00 

+3100.0 

+00 

+6918.0 

+ 00 



TABLE I - SAMPLE PROBLEM INPUT - Concluded 


+ 3pnu *i) 
+ 3t>ou • 

+ ^rtOtj « it 

+41MU.M 
+ uh no • o 


(b) Fortran data card listing 


+u n +7uSu • u 
+UH+7UPU . 0 
+ u f ' + r* 7 u • o 
+uo+Rbnn . n 

+u r ‘+q l ,nu • u 


+uo+33no. o 
+on+3bno.u 

+ 011 + V+ntj . u 
+uo+42on. o 

+on+4bno.o 


+00+7175.0 

+on+7b30.o 

+on+8i?o.n 

+nn+R700 . n 
+00+9150.0 


+QO+34O0 . 0 
+00+3700.0 

+on+anno . o 

+00+4300.0 

+on+4h0o.o 


+00+7350. 0 
+00+7800.0 
+00+8300.0 
+00+885 0.0 
+00+9270.00 


+ 00 
+ 00 
+ 00 
+ 00 
+ 00 



TABLE II.- SAMPLE PROBLEM OUTPUT 


TYPICAL CHARRING ARL ATDR - TEST CASE - 4/6/65 DONALD M. CURRY 


INPUT DATA. 


TIME LIMIT = 6. OOUO E 02 INITIAL TlME=0. 


NPTT= 2 


TIME 

0 . 

6.00-J0E 02 


TIME STFP 
I .OOUOE-OI 
I . 3CO0E-01 


PRINT CONTROL 
i DO 
100 


F CON V= l.OODPQE 00 


FRAD- UOOOOOE J 0 


TYPICAL CHARRING ABLATION MATERIAL PROPERTIES 


T ABL = 1 . 0 6 0 0 j E J 3 
F BLOW* 0. 

H V- 2. 50 COO C 02 
C HAR C = 4. 30 00 DC-01 

NP- 31 NKC= 


TCHAR= 1.46COOE 03 
EMV = 6. 5C000E-01 
VP T = 0. 

ABL K = 7.00000E-02 
2 N C P C = 2 Nl 


TR E C = 0. 

EMC= 7 . 500C0E-0 l 

f v* l.nooooe oo 

AB LC = 4 . 30000E-0 1 
9 NC PV = 2 


RHOV- 3.40000E 01 
H300= 1 . 29060 E 02 
I V= 0 . 


NR EC= 2 


RHOC= 2.00000E 01 
V L = 1.50000E 00 
CH ARK= 1 . 2 OOOOE- 0 1 


TEMPERATURE 
3.6000OF 02 
4. 60000 E 02 
5 . 60000 E 02 
6.60000E 0? 
7 . 60000 E 02 
8 . 60000 E 0? 
9 • 6 0 0 0 0 E 02 
1 ■ 06000 E 03 
1 . 16 300 E 03 


THERMAL 
CONDUCTIVITY 
b. 8COO0E-O2 
6*5000 OE- 02 
6 . 5 5 0 0 J E - U 2 
6.6C00DE-02 
6. 72000E-02 
6.84000E-02 
6.90000E-C2 
7 • Uu 000 t -02 
7.D000JE-02 


VIRGIN MATERI 


TEMPERATURE 
3.60000E 02 
1 • 1 OOOOE 03 


SPECIFIC 

HEAT 

H.30000E-01 
4 .30000E-01 


TEMPERATURE 
1 . 46000 E 0 3 
l.OOOOOC 04 


THERMAL 
CONDUCT I VI I Y 
1 . 20 00 OE- 01 
1 .2000JE-01 


CHAR MATERIAL 


TEMPERATURF 
1.46D00E 03 
1 . C0GOGE 04 


SPEC I F IC 
HEAT 

t.30U0’Jt vn 
4. 30000E-01 


SURFACE RECESSION TABLE 


TIMF 

0 . 

6 * OOOOOE 02 


SK - IN/SEC 
9 . OOQOOE-04 
9.0C000E-04 


NO TRAJECTORY - 0=95 BTU/SEC-SQFT 


NO. 


2 


OF TRAJECTORY POINTS 



TABLE II.- SAMPLE PROBLEM OUTPUT - Continued 


TIME 

C. 

6.OOC1GOE 02 


0 CONVECTIVE 
*.DUuUUfc U 1 

9.50000E 01 


0 RADIATIVE 

0 . 

0 . 


VELOCITY 
2.9260UC uA 
2.92500E OA 


PROPERTIES OF BACKUP STRUCTURE 


NO. OF MATERIALS IN BACK-UP SH IELD= 1 
TOTAL DUMBER OF NUDES IN BACK-UP SH I ELD= 3 
THICKNESS OF LACK-UP 3HIEL0= i.cooooe-oi 


BACKUP 

MAT 

ER I A L 

0.1 INCHES THICK 




T II F f* M ft 1 


SPECIFIC 

TEMPERATURE 

CUN DUG 1 1 V 1 1 r 

TEMPERATURE 

nEAT 

3. 6 00 O'j E 

0? 

6. BOO 00 E- 02 

3.60C00E 02 

A.30000E — 01 

A • 60000 F 

0? 

6 . 5 C 0 0 0 E - 0 2 

1.10 CODE 03 

A . 300C0E-0 1 

5. 60000 E 

02 

6. *500 ) E-02 



6.6000Tb 

02 

6.6C0UUE-02 



7. 60000 E 

02 

6 . 7 ?O00 c — 02 



8 . 6 0 0 0 ' 1 F 

02 

6.dAOOJ*>02 



9. 60000 E 

02 

6. JCOOjE-O? 



1.060 ODt 

03 

7. OOOOOfc-02 



1 . 160U0C 

03 

f . ''Q n 00Z~02 




MATERIAL 

DENSITY 

THICKNESS 

EMI SSI VI TY 
FRONT BACK 

NGDES/MATERIAL 

1 

3.AOOOE 01 

i. 00 DOE-01 

9 . 0000c-0 1 9 . 0000 E-0 1 

3.0000E 00 

ADDITIONAL DATA 

FUR IMDIVIDUAL MAlERlALb 

IN BACKUP STRUCTURE 


MATERIAL 

I 

FILM COEFFICIENT 

n. 

GAP THICKNESS FTEST 

0. 0. 

BTEST 

0. 


HEAT TkANSFLR Til CABIN ENVIRONMENT - HENV'O.O 
TtMPtKATUKI> 6 • 600 00 E 02 FILM COEFFICIFNl= 0. VIEW FACTOR= 0. 


INITIAL TEMPERATURE IS CONSTANT 

TEMPFRATIihf DISTRIBUTION IN HEAT SHTRO IS UNIFORM AND FOUAL TO 5.3000E 02 


Q LO S T= 0. 


99 



TABLE II.- SAMILE PROBLEM OUTPUT - Concluded 


OUTPUT U AT A . 


TIME= 9.90000E 00 QCUNVECTIVE= 9.50000E 01 QRAOIATIV£= 0. VELOCITY= 2.92500E 04 

GAS ABLATION RATE= 0. CHAR ABLATION R ATE = 5.40000E 00 TOTAL ABLATION R ATE= 5.40000E 30 

RECESSION DEPTHS 9.00000E-03 QHOT WALL= 8.99202E 01 

TEMPERATURE DISTRIBUTION IN HEAT SHIELD AT THE END OF THE TIME STEP, T= l.OODOOE 01 SECONDS 


TEMPERATURE 

DISTRIBUTION IN 

THE 

ABLATING MATERIAL 







3.79022E 

03 

2.33046E 

03 

1.06300E 

03 

6.38075E 

02 

5.47600E 

02 

5.32434E 

02 

5.30239b 

02 

5.30030E 

02 

5 • 30002 E 

02 

5.30000E 

02 

5.30000E 

02 

5.30000E 

02 

5.30000E 

02 

5.30000E 

02 

5.29999E 

02 

5.30000E 

02 

5.30000E 

02 

5. 3300 OE 

02 

5.30000E 

02 

5.30000E 

02 

5.30000E 

02 

5.30000E 

02 

5.30000E 

02 

5.30000E 

02 

5.30000E 
5 .299996 

02 

02 

5.30000E 

02 

5.30000E 

02 

5.29999E 

02 

5.30000E 

02 

5.300006 

02 


TEMPERATURE 

DISTRIBUTION IN 

THE 

BACK-UP STRUCTURE 

5.29999E 

02 5.30000E 

02 

5.30000E 02 


T I ME= 1.99000E 01 QCONVECTIVE= 9.50000E OL QRADI AT I VE= 0. VELOCITY= 2.92500F 04 

GAS ABLATION R A T b= 2.90720b 01 CHAR ABLATION RATE= 5.40000E 00 TOTAL ABLATION RA TE= 3.44720E 01 

RECESSION DEPTH= 1.80000E-02 QHOT WALL= 8.99173E 01 

TEMPERATURE DISTRIBUTION IN HEAT SHIELO AT THE END OF THE TIME STEP, T= 2.00000E Oi SECONDS 


TEMPERATURE 

DISTRIBUTION IN 

THE 

ABLATING MATERIAL 







3.80935E 

03 

2.8767BE 

03 

1.64575E 03 

9.77202E 

02 

6.72698E 

02 

5.66204E 

02 

5.37857E 

02 

5.31496E 

02 

5.30254E 02 

5.30038E 

02 

5.30005E 

02 

5.33300E 

02 

5.29999E 

02 

5.29999E 

02 

5.29999E 02 

5.29999E 

02 

5.29999E 

02 

5.29999E 

02 

5 .299996 

02 

5.29999E 

02 

5.29999E 02 

5.29999E 

02 

5.29999E 

0? 

5.29999E 

02 

5.29999E 

02 

5.29999E 

02 

5.29999E 02 

5.29999E 

02 

5.29999E 

02 

5.29999E 

02 

5.29999E 

02 










TEMPERATURE 

DISTRIBUTION IN 

THE 

BACK-UP STRUCTURE 







5.29999E 

02 

5.29999E 

02 

5.29999E 02 







TIME = 2.99000E 

01 

OCONVECTI VE= 

9. 

50000E 01 QR ADI AT I VE= 0. 


VELOC I TY= 

2.92500E 04 


GAS ABLATION RATE= 

1.25330E 01 

CHAR ABLATION R A TE = 5 

• 40 00 OE 00 

TOTAL 

ABLATION 

RA TF = 

1 . 79 33 OE 

31 

RECESSION DEPTH 

= 2 

.700006-02 QHOT 

WALL= 8.98427E 01 








100 













Temperature, 




Exact 

Locatic 

solution 

in. 

O 

0.2 

□ 

0.4 

O 

0.6 

A 

0.8 

STAB 





Lot 


3200 1 - 




O Test 
point 

v 


Thermocouple 
depth, in. 

0.1 
0.3 
0.5 
0. 75 
1.0 

STAB II results 


charring ablator. 



Surface recession, in. 



Time, sec 


Figure 8. - Plot program surface recession curve from typical charring 
ablator test case. 





Time, sec 


Figure 9. - Plot program bondline temperature curve from 
typical charring ablator test case. 



Distance, in 



Time, sec 


Figure 10. -Plot program 1060°R and 1460°R isotherm curves from 
typical charring ablator test case. 
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I 


I Mill HI II l■ll■ll■■llll Min III 


"The aeronautical and space activities of the United States shall be 
conducted so as to contribute ... to the expansion of human knowl- 
edge of phenomena in the atmosphere and space. The Administration 
shall provide for the widest practicable and appropriate dissemination 
of information concerning its activities and the results thereof.” 

— National Aeronautics and Space Act of 1958 


NASA SCIENTIFIC AND TECHNICAL PUBLICATIONS 


TECHNICAL REPORTS: Scientific and technical information considered 

important, complete, and a lasting contribution to existing knowledge. 

TECHNICAL NOTES: Information less broad in scope but nevertheless 

of importance as a contribution to existing knowledge. 

TECHNICAL MEMORANDUMS: Information receiving limited distri- 

bution because of preliminary data, security classification, or other reasons. 

CONTRACTOR REPORTS: Technical information generated in con- 

nection with a NASA contract or grant and released under NASA auspices. 

TECHNICAL TRANSLATIONS: Information published in a foreign 

language considered to merit NASA distribution in English. 

TECHNICAL REPRINTS: Information derived from NASA activities 

and initially published in the form of journal articles. 

SPECIAL PUBLICATIONS: Information derived from or of value to 

NASA activities but not necessarily reporting the results of individual 
NASA-programmed scientific efforts. Publications include conference 
proceedings, monographs, data compilations, handbooks, sourcebooks, 
and special bibliographies. 


Details on the availability of these publications may be obtained from: 
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NATIONAL AERONAUTICS AND SPACE ADMINISTRATION 

Washington, D.C. 20546 



